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The  Nonlinear  Control  Theory  of  Complex 
Mechanical  Systems 

Grant  AFOSR  90-0226 


FINAL  REPORT 


Abstract 

This  report  provides  a  technical  overview  of  research  on  the  nonlinear  control  the¬ 
ory  of  complex  mechanical  systems,  conducted  under  support  of  the  U.S.  Air  Force 
Office  of  Scientific  Research.  The  period  of  support  has  been  June  15,  1990  through 
December  14,  1995,  The  research  has  had  both  theoretical  and  experimental  compo¬ 
nents  aimed  at  development  of  new  approaches  for  controlling  the  motions  of  complex 
mechanical  systems.  The  work  has  been  principally  concerned  with  classes  of  physical 
systems  wherein  the  salient  dynamical  features  cannot  be  understood  in  terms  of  linear 
models.  Applications  of  interest  include  the  control  of  molecular  dynamics,  microelec¬ 
tromechanisms,  aerospace  structures,  rotating  shafts,  turbine  dynamics,  etc.  Common 
elements  in  models  of  such  systems  form  the  basis  of  a  general  theory  of  control. 

The  principal  focus  of  this  report  is  on  control  strategies  using  tuned  oscillatory 
forcing.  Two  papers  in  particular — “Stable  Average  Motions  of  Mechanical  Systems 
Subject  to  Periodic  Forcing”  and  “Energy  Methods  for  Stability  of  Bilinear  Systems 
with  Oscillatory  Inputs”  provide  a  detailed  account  of  some  of  our  research.  The 
broader  scope  of  our  research  is  indicated  by  the  list  of  references  which  concludes  the 
report. 
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1  Introduction  and  Report  Summary 

During  the  period  of  performance  covered  by  this  report,  basic  research  concerned  with  the 
nonlinear  control  theory  of  complex  mechanical  systems  was  carried  out  at  Boston  Univer¬ 
sity.  Results  of  the  research  have  been  reported  in  (i)  28  technical  papers  appearing  in 
the  open  literature  (references  [15-42]  cited  at  the  end  of  the  report)  and  in  (ii)  numerous 
technical  presentations  at  conferences,  AFOSR  contractors’  review  meetings,  and  lectures  at 
educational  institutions  throughout  the  world.  In  May,  1993,  two  doctoral  students  who  had 
been  supported  on  the  subject  grant  received  Ph.D.  degrees  from  Boston  University.  These 
were: 

1.  Daniel  P.  Martin,  B.U.  Ph.D.  Thesis:  “Mathematical  Methods  for  Problems  of  Kine¬ 
matic  Redundancy  in  Robotics,”  and 

2.  Danbing  Seto,  B.U.  Ph.D.  Thesis:  “Stabilization  Problems  in  the  Control  of  Super- 
Articulated  Mechanical  Systems.” 

The  research  covered  by  this  report  has  been  focused  on  several  different  problem  areas: 

(a)  Robotic  devices  with  redundant  degrees  of  freedom  (Refs.  [16],  [18], [22], [27-29], [32]). 

Some  of  the  early  research  sponsored  under  the  subject  grant  dealt  with  planning 
motions  for  robotic  mechanisms  which  simultaneously  exhibit  features  of  joint  elasticity 
and  kinematic  redundancy  ([27],  [28],  and  [29]).  It  has  been  shown  that  for  idealized 
models  of  mechanisms  with  certain  congenial  geometries,  it  is  possible  to  plan  motions 
which  never  store  energy  in  their  elastic  degrees  of  freedom.  Such  motions  follow 
acceleration  constraints  which  are  most  naturally  described  in  a  coordinate  system 
attached  to  some  component  of  the  mechanism.  Results  to  date  point  to  two  lines  for 
further  study:  (i)  It  is  apparently  only  highly  idealized  (non-physical)  models  which 
admit  the  possibility  of  moving  without  storing  energy.  We  wish  to  develop  a  rational 
basis  for  associating  such  models  to  actual  systems  with  the  aim  of  ascertaining  whether 
“ideal”  motions  which  result  in  no  energy  storage  in  the  model  also  prescribe  motions 
which  minimize  energy  storage  in  the  actual  system,  (ii)  Current  research  by  us  and 
others  seeks  to  develop  control  strategies  which  implement  the  motion  plans  that  have 
been  proposed. 

(b)  Modeling  and  control  of  the  rotational  dynamics  of  complex  multibody  systems  and 

mixed  structures  (Refs.  [20], [24], [33]). 

Our  early  work  in  developing  a  control  theory  of  mechanical  systems  concentrated  on 
development  of  models  of  rotating  mixed  structures  having  multiple  interconnected 
elastic  and  rigid  components.  Large  mechanically  complex  earth  satellites  were  the 
primary  applications  focus  of  this  research.  While  there  had  been  a  great  deal  of 
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prior  work  on  the  control  of  space  structures,  there  was  little  literature  treating  the 
application  of  ideas  from  geometric  mechanics  and  nonlinear  control  theory  to  problems 
where  such  ideas  were  of  apparent  importance.  In  [3],  we  proposed  a  general  Lagrangian 
formulation  for  the  dynamics  complex  mechanisms  capturing  the  general  effects  of 
inertial  forces  created  by  spatial  rotations.  Using  Rayleigh  type  dissipation  models, 
we  developed  a  simple  geometric  formalism  for  modeling  rotating  viscoelastic  systems. 
This  theory  was  applied  to  analyze  the  dynamic  behavior  of  a  simple  spatially  rotating 
structure  consisting  of  a  rigid  body  with  a  viscoelastic  beam  attachment.  Asymptotic 
steady  state  behavior  of  the  model  was  determined,  and  a  complete  stability  and 
bifurcation  analysis  was  given  for  a  closely  related  model  having  only  one  rotational 
degree  of  freedom.  It  is  noted  that  parallel  research  on  the  modeling  and  analysis  of 
mixed  structure  dynamics  was  also  carried  out  with  greater  emphasis  on  Hamiltonian 
methods  by  Krishnaprasad,  Marsden  and  others. 

(c)  Adaptive  and  robust  control  of  super-articulated  (underactuated)  mechanical  systems 

(Refs.  [31],[34],[36]). 

Some  of  the  research  cited  in  the  references  at  the  end  of  the  report  has  been  aimed 
at  classifying  multibody  systems  in  terms  of  their  topology  and  mechanical  intercon¬ 
nections  (joint  types  between  bodies).  (See,  e.g.,  [31]  and  [36].)  By  assigning  a  certain 
digraph,  which  we  call  a  control  flow  diagram  (CFD),  to  each  super-articulated  mechan¬ 
ical  system,  we  develop  a  formal  notion  of  control  complexity  which  provides  informa¬ 
tion  on  how  control  design  may  be  usefully  approached.  Certain  “chain”  systems,  for 
instance,  are  shown  to  be  feedback  controllable  using  a  backstepping  controller  design 
whose  stability  properties  are  easily  deduced  from  the  backstepping  algorithm  itself. 
More  complex  structures  call  for  greater  design  effort,  but  the  graphical  methods  that 
have  been  studied  are  useful  in  prescribing  decoupling  control  laws  which  effectively 
transform  certain  structures  into  systems  of  chains  to  which  our  design  methods  apply 
directly.  The  CFD  also  points  to  a  classification  of  mechanical  systems,  for  which  we 
may  design  stable  adaptive  feedback  control  laws  to  treat  parametric  uncertainties. 
(See  [34].) 

(d)  Stability  and  control  of  mechanical  systems  using  oscillatory  actuator  inputs  (Refs. 

[17],[25],[26],[30],[37-41]). 

Motivated  by  our  early  work  on  the  control  and  bifurcation  theory  of  equilibrium  rota¬ 
tions  in  multibody  systems  (See  [1-6],  [11],  and  [24].),  we  have  turned  our  attention  to 
the  broader  issue  of  open-loop  control  designs  tailored  to  exploiting  key  the  intrinsic 
dynamical  behavior  in  systems  of  interest.  More  specifically,  we  have  been  studying  the 
use  of  oscillatory  forcing  to  produce  stable  motions  in  a  variety  of  mechanical  systems, 
including  rotating  and  vibrating  kinematic  chains  as  well  as  in  the  axial  compressor 
problems  discussed  below.  The  mathematical  mechanism  underlying  this  approach  to 
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stabilization  involves  the  nonholonomic  relationship  between  input  and  state  variables 
in  a  fundamental  way.  (See  [40]  and  [37]  for  details.)  Using  various  extensions  of 
classical  averaging  methods,  it  has  been  shown  that  the  stabilizing  effects  of  oscillatory 
forcing  may  be  identified  and  analyzed  in  terms  of  an  energy-like  quantity  call  the 
averaged  potential.  (See  [30]  [40],  and  [37].)  Under  certain  reasonable  assumptions, 
forcing  by  periodic  inputs  produces  dynamic  responses  in  which  states  of  the  system 
are  confined  to  neighborhoods  of  local  minima  of  the  averaged  potential.  While  use  of 
the  averaged  potential  to  analyze  the  stabilizing  effect  that  can  be  produced  by  appro¬ 
priately  designed  oscillatory  forcing  of  conservative  (Hamiltonian)  mechanical  systems 
highlights  the  geometric  character  of  the  phenomenon,  we  have  also  begun  looking  at 
this  approach  for  influencing  the  stability  characteristics  of  dissipative  systems.  The 
main  body  of  this  report  consists  of  two  papers  we  have  published  on  the  use  of  oscil¬ 
latory  forcing  for  stabilizing  finite  dimensional  mechanical  systems. 
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2  Stable  Average  Motion  of  Mechanical  Systems  Sub¬ 
ject  to  Periodic  Forcing 
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Stable  Average  Motions  of  Mechanical  Systems  Subject 

to  Periodic  Forcing 

J.  Baillieul  * 


Abstract:  There  is  a  rapidly  growing  body  of  literature  devoted  to  the  study  of  anholonomy  in  the 
motions  of  nonlinear  systems  forced  by  periodic  inputs.  While  a  great  deal  of  recent  work  has  treated 
the  case  of  systems  without  ^‘drift”,  the  emphasis  in  our  treatment  will  be  the  controlled  dynamics  of 
mechanical  systems  in  which  “drift”  terms  play  an  important  role.  Our  main  application  will  be  to 
controlling  the  dynamics  of  “super-articulated”  mechanical  systems — systems  having  fewer  controls 
than  configuration-space  dimensions.  Using  classical  averaging  theory,  it  is  shown  that  the  stability 
of  motion  for  such  forced  systems  may  be  analyzed  using  energy  methods  together  with  the  adroit 
introduction  of  dissipation  into  the  models.  For  systems  where  the  control  inputs  are  applied  directly 
only  to  cyclic  coordinates,  we  define  a  simple  but  important  quantity  called  the  averaged  potential. 
It  is  shown  that  under  certain  reasonable  assumptions,  forcing  by  periodic  inputs  produces  dynamic 
responses  which  are  confined  to  neighborhoods  in  the  phase  space  associated  with  local  minima  of 
the  averaged  potential.  We  conclude  with  a  discussion  of  linearization  of  the  forced  system  about 
critical  points  of  the  averaged  potential.  While  the  critical  points  of  the  averaged  potential  may 
not  correspond  to  equilibrium  points  of  the  original  forced  system,  averaging  theory  implies  that 
trajectories  of  the  forced  system  will  “hover”  around  such  points.  For  an  example  problem  in 
which  there  is  an  energy-like  quantity  which  is  conserved,  it  is  shown  by  purely  geometric  means 
that  Lyapunov  stability  may  be  deduced  from  a  critical  point  analysis  of  the  averaged  potential. 
Hence,  although  the  role  of  the  averaged  potential  in  our  stability  analysis  is  initially  established 
by  a  classical  averaging  argument  applied  to  an  appropriate  dissipative  system,  it  is  shown  to  be  a 
geometric  quantity  which  describes  the  effects  of  anholonomy  and  is  independent  of  stability  in  the 
presence  of  dissipation. 


1.  Introduction 


Consider  a  cart  which  is  free  to  move  on  a  frictionless  track  and  to  which  there  is  attached  a 
pendulum  with  a  frictionless  hinge  whose  axis  is  aligned  with  gravity.  (Thus  we  have  a  two  degree 
of  freedom  mechanical  system  in  which  there  are  no  potential  forces  which  influence  the  motion.) 

*The  author  gratefully  acknowledges  the  support  of  the  U.  S.  Air  Force  through  grant  AFOSR-90-0226. 
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Figure  1:  Phase  portrait  of  (2)  under  periodic  forcing  x  =  sin(9t). 

Suppose  the  mass  of  the  cart  is  and  the  pendulum  consists  of  a  rigid  massless  rod  of  length  f 
together  with  a  tip  of  mass  rrip.  Assume  further  that  we  may  apply  a  control  force  ii  to  move  the 
cart  along  its  track,  and  no  other  exogenous  forces  are  applied.  Let  x  denote  the  distance  (from 
some  arbitrary  reference  point)  that  the  cart  has  moved  along  the  track,  and  let  0  denote  the  angle 
of  the  pendulum  with  respect  to  the  positive  direction  of  the  track,  x  and  0  together  comprise  a 
S3^stem  of  generalized  coordinates  in  terms  of  which  we  may  write  the  dj^namics  of  the  sj^^stem 

(7Tic  +  ?np)x  —  mpCO  sin  0  —  nipi  cos  69^  =  u  (1) 

£6  —  X  sin  6  =  0.  (2) 

With  the  single  scalar  control  entering  the  dynamics  in  this  way,  the  cart-pendulum  system  is 
a  member  of  the  class  of  mechanical  systems  which  we  have  called  super-articulated.  (A  super- 
articulated  mechanical  systems  is  a  controlled  mechanical  system  wherein  the  number  of  controls 
is  strictly  less  than  the  dimension  of  the  configuration  space.)  Broadly  speaking,  the  goal  of  the 
present  paper  is  to  develop  the  elements  of  a  control  theory  of  super-articulated  mechanical  systems 
which  will  explain  the  behavior  of  (l)-(2)  which  is  illustrated  in  Figure  1.  This  shows  a  phase 
portrait  of  (2)  under  periodic  forcing  x  =  — L  cosat  with  initial  conditions  {9,9)  =  (•,  •).  The  point 
(^9, 9)  =  (0, 0)  appears  to  be  Lyapunov  stable,  and  a  proof  of  the  fact  that  it  is  will  emerge  from  our 
theory. 

The  remainder  of  the  paper  is  organized  as  follows.  In  Section  2,  we  define  what  is  meant  by  a 
generalized  mechanical  control  system.  This  includes  both  classical  finite-dimensional  mechanical 
control  systems  (e.g.  robots  with  rigid  components)  and  “super-articulated”  mechanical  systems 
(e.g.  robots  with  elastic  joints)  as  special  cases.  Various  examples  are  discussed.  In  Section  3, 
we  apply  classical  averaging  theory  to  the  study  of  periodically  forced  generalized  control  systems. 
The  averaged  potential  is  defined  in  Section  4.  This  function  is  shown  to  be  useful  in  characterizing 
the  stability  of  motions  in  which  classical  potential  forces  interact  with  inertial  forces  created  by 
the  motion  produced  by  the  periodic  forcing.  While  the  properties  of  the  averaged  potential  are 
illuminated  using  classical  averaging  theory  in  Section  4,  Section  5  pursues  a  geometric  analysis 
which  shows  that  the  averaged  potential  provides  a  simple  measure  of  nonholonomic  effects  produced 
by  periodic  forcing  of  the  class  of  systems  under  study. 


STABLE  AVERAGE  MOTIONS  OF  MECHANICAL  SYSTEMS 
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This  paper  presents  a  geometrical  framework  within  which  the  types  of  mechanics  problems  typically 
treated  by  classical  averaging  theory  may  also  be  solved  using  methods  suggested  by  geometric  non¬ 
linear  control  theory.  Advantages  of  this  approach  include  a  sharper  characterization  of  dynamical 
behavior  and  more  direct  treatment  of  conservative  (Hamiltonian)  systems.  Our  geometric  analysis 
in  Section  5  is  aimed  at  understanding  the  effects  of  anholonomy  in  certain  periodically  forced  non¬ 
linear  systems,  and  in  this  work  we  have  drawn  inspiration  from  the  recent  work  of  Sussmann  (e.g. 
[17],  [18]),  Sastry  ([13]),  their  coworkers,  and  others  ([10]).  The  problems  we  treat,  however,  differ 
significantly  from  “kinematic  control”  problems  in  that  “drift”  terms  are  present  in  the  dynamics 
we  treat,  rendering  both  our  methods  and  results  quite  different  from  those  reported  in  [17],  [18], 
and  [13]. 


2.  Generalized  Mechanical  Control  Systems 


Following  Fliess  ([7])  we  state  the  following: 


Definition  1  A  generalized  second-order  control  system  is  prescribed  by  an  evolution  equation  of 
the  form 

q{t)  =  F{q{t),  q{t);  u,u,...,  u^"’) 

■where  q  takes  values  in  a  configuration  manifold  M  and  control-input  curves  u[t)  are  bounded  and 
piecewise  analytic  with  values  in  a  control-input  manifold  U . 


NOTE:  Assuming  that  the  control  u  is  onl}'  piecewise  analytic  means  that  the  derivatives  ii,. . .  may 
not  be  defined  at  all  times  t.  This  is  a  technicality  which  may  be  addressed  by  arbitrarily  redefining 
u{t)  and  it  derivatives  at  points  of  ambiguity. 

Since  we  shall  be  interested  in  local  behavior,  there  is  no  loss  of  generality  in  assuming  M  and  U 
are  finite  dimensional  vector  spaces  of  appropriate  dimensions.  As  in  [7],  we  further  assume  that 
for  each  (9,  q)  €  TM^  F{q^  q]u,. . is  a  polynomial  in  u  and  its  derivatives. 


Definition  2  A  generalized  mechanical  control  system  is  a  generalized  second  order  control  system 
which  has  second-order  polynomial  dependence  jointly  on  x  and  u  and  which  has  polynomial  depen¬ 
dence  on  ii  which  is  at  most  first  order.  Hence  the  equations  of  a  mechanical  control  system  take 
the  following  form: 

Q  =  9o{q,  q)  +  gi{q)u  -f  g2{q,  q)u  +  g3{q)u^'^^.  (3) 

where  denotes  the  symmetric  second  tensor  power  of  the  vector  u  (See  [6]  for  details!),  and 
where  for  each  q,  go{q,q)  has  at  most  second  order  polynomial  dependence  on  q,  and  g2{q,q)  has  at 
most  first  order  polynomial  dependence  on  q. 
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Remark  1 
as 


To  emphasize  the  polynomial  dependence  on  q  we  write  out  the  velocity  terms  in  (3) 

9o{q,  q)  =  gM  +  9o{g)g  +  gliqW^K  (4) 

92{q,q)  ^  9M  + 9l{q)q-  (5) 


It  is  useful  to  examine  some  examples  of  mechanical  control  systems. 


Example  1  (Cart-pendulum  1)  Equation  (2)  may  be  viewed  as  a  specialization  of  (3)  in  which  0 
is  the  configuration  variable,  x  is  the  control  (i.e.  u  =  i),  go  =  g\  —  go  =  0,  and  g2  is  a  function  of 
0  alone  with  g2{0)  =  sm6. 


Example  2  (Cart-pendulum  2)  Equation  (2)  may  also  be  viewed  as  a  specialization  of  (3)  in  which 
9  is  the  configuration  variable,  x  the  control  (i.e.  u  =  x),  go  =  g2  =  go  =  0)  ^.nd  gi{9)  =  —  sin0. 


Remark  2  The  structure  (3)  is  sufficiently  rich  that  there  may  be  multiple  ways  in  which  to 
represent  a  particular  physical  system  of  interest.  The  choice  of  representation  will  be  guided  by 
the  advantages  afforded  in  applying  the  theory. 


Example  3  (Lagrangian  control  systems)  Let  M  be  a  real  analytic  differentiable  manifold  of  di¬ 
mension  n.  Following  Nijmeijer  and  van  der  Schaft  ([14]),  we  define  a  Lagrangian  control  system  on 
M  to  be  a  dynamical  system  with  inputs  whose  equations  of  motion  are  prescribed  by  applying  the 
Euler- Lagrange  operator  to  a  function  L  :  TM  x  U  ^  R,  L  =  L{q,q-,  u),  whose  dependence  on  the 
configuration  q,  the  velocity  q,  and  the  control  input  u  is  smooth.  The  set  U  in  which  the  control 
input  functions  u(-)  take  values  is  assumed  to  be  a  closed  convex  subset  of  R™  which  we  further 
assume  is  symmetric  with  respect  to  the  origin.  (I.e.,  u  G  17  —u  G  U.) 


A  nondegenerate  Lagrangian  is  one  for  which  nonsingular  for  all  {q,  q\  u)  G  TM  xU. 

The  equations  of  motion  for  a  control  system  arising  from  a  nondegenerate  Lagrangian  may  be 
written 

.  dH  .  dL. 


dq^ 


dqdq  dudq  dq 


Defined  in  this  way,  a  Lagrangian  control  system  is  a  generalized  mechanical  control  system  precisely 
when  L  has  the  proper  polynomial  dependence  on  the  variable  u. 


A  wide  class  of  interesting  mechanical  control  systems  arise  from  Lagrangian  functions  of  the  form 
7'(9,4;  u)  =  \q^M{q)q  -  V{q)  +  q'^Bu  where  q  G  R”,  u  G  R™,  and  R  is  an  n  x  m  matrix. 
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Example  4  (Reduced-order  symmetric  Lagrangian  systems)  Consider  a  Lagrangian  control  system 
with  configuration  variables  (91,92)  €  R"‘  x  R"L  The  variable  91  will  be  called  cyclic  if  it  does 
not  enter  into  the  Lagrangian  when  u  =  0,  i.e.  if  ^(91, 92, 9i,  92;  0)  =  0.  We  shall  be  interested  in 
Lagrangian  control  systems  with  cyclic  variables  in  which  it  is  precisely  the  cyclic  variables  which 
may  be  directly  controlled.  Specifically,  we  shall  consider  systems  of  the  form 


^(9i>?2,9i,92W)  =  ^(9m9D 


(  "^(92) 
\  ^(92) 


(6) 


where  dim 91  =  ni,  dim 92  =  ?^2,  ^rid  m  and  M  are  symmetric  positive  definite  matrices  of  dimension 
72]  X  Ui  and  722  x  722  respectively.  To  emphasize  the  distinguished  role  to  be  played  by  the  cyclic 
variable,  we  write  v  =  qi. 


Applying  the  usual  Euler-Lagrange  operator  to  this  function,  we  find  the  equations  of  motion  are 


^(722(92)2;  +  ^(92)^92)  =  72 

92  +  M(92)~^(r(92,92)  +  ^  +  ^(92)22  +  ^(92,92)2^  -  ^22^1^22)  =  0 


dq2 


dq2 


where  the  A;-th  component  (1  <  A:  <  722) 


1  ^,dMki  ,  dMkj  dMii 


r(92,92)fc  -  X  ^(-  j  O 

2^  dqi  dq\  C>92 


F  )^292i 


and  the  A:-th  component 


(7) 

(8) 


Noting  that  u  and  v  have  the  same  dimensions  and  that  722(9)  2s  nonsingular  for  all  9,  it  is  an 
immediate  consequence  of  existence  and  uniqueness  theorems  for  ordinary  differential  equations 
that  for  any  continuous  trajectory  92(t),  any  smooth  trajectory  v{-)  may  be  generated  by  (7)  by 
means  of  an  appropriate  choice  of  control  input  u{-).  This  means  that  in  principle,  it  is  immaterial 
in  controlling  (7)-(8)  whether  we  view  u  as  the  control  with  v  being  determined  by  integrating 
the  differential  equation  (7)  or  take  v  to  be  the  control  with  u  then  simply  determined  by  (7).  If 
V  is  viewed  as  the  control  input,  (8)  clearly  has  the  form  of  a  mechanical  control  system  (3).  It 
is  interesting  to  note  that  (8)  itself  prescribes  the  dynamics  of  a  Lagrangian  control  system  with 
associated  reduced  Lagrangian 

L{q2, 92;  27)  =  ^9^ -^(92)92  +  9j^(92)27  -  K(9;  22), 

where  K  is  the  augmented  potential  defined  by  Va{q',v)  =  V{q)  —  ^v^m{q2)v.  The  terms  ^(92,92)22 
which  are  linear  in  the  generalized  velocities  are  called  gyroscopic  terms.  Many  physical  systems 
admit  this  type  of  reduction,  including  the  cart-pendulum  system  (l)-(2)  and  the  rotating  kinematic 
chain  systems  studied  in  [2]  and  [12].  For  more  details  on  the  role  played  by  gyroscopic  terms  in 
the  dynamics  of  mechanical  systems,  the  reader  is  referred  to  [3]  and  [9]. 
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3.  Averaging  and  the  Stability  of  Motion 


Although  many  interesting  physical  systems  may  be  described  by  models  of  the  form  (3),  the 
control  theory  of  such  systems  is  presently  incomplete.  While  it  is  generally  unclear  where  to  turn 
for  inspiration  in  classical  geometric  nonlinear  control  theory,  recent  work  by  Murray  and  Sastry, 
[13],  suggests  that  the  geometry  of  noncommuting  vector-fields  may  be  systematically  exploited  by 
using  carefully  selected  time-periodic  inputs  to  induce  desired  motions.  (See  also. related  work  by 
Sussmann,  [18].)  Such  results  are  not  immediately  applicable  to  systems  of  the  form  (3),  however, 
and  it  is  our  goal  in  the  remainder  of  the  paper  to  develop  the  appropriate  theory.  In  the  present 
section,  our  principal  mathematical  tool  will  be  classical  averaging  theory. 

Suppose  the  input  u(-)  appearing  in  (3)  is  a  periodic  function  of  time  with  fundamental  period 

T  >  0.  To  put  (3)  into  a  form  to  which  standard  results  in  averaging  theory  may  be  applied,  we 

convert  to  first  order  form  and  rewrite  the  equations  in  the  time  scale  r  =  tjT.  (Cf.  [8],  p.  166.) 
Letting  ?/(r)  =  x{t),  the  resulting  differential  equations  take  the  form 

—  ( \  =  T  ■  (  ri')  ,  (9) 

dr  Vy2/  \go{yi,y2)  +  9i{yi)'ds  +  g2{yi,y2)us  +  gsiviWs  /  ' 

where  Us{r)  =  u{Tt),  Us(t)  =  u{Tt).  As  the  frequency  of  the  forcing  is  varied,  we  find  the 
amplitude  of  ii,  grows  like  (Think  of  varying  the  frequency  of  a  sinusoid.)  Hence,  averaging 
theory  does  not  immediately  apply  unless  gi  =  0.  Fortunately,  a  simple  time- varying  coordinate 
transformation  eliminates  the  dependence  on  the  derivative  of  the  input. 


Definition  3  Given  the  mechanical  system  (3),  the  L-transformation  is  defined  by  the  equation 


X2J  \q-gi{q)uJ 


(10) 


Important  features  of  this  transformation  are  given  in  the  following. 


Proposition  1  (i)  The  L-transformation  is  invej'tible.  (ii)  In  [x\^X2)-coordinates  (3)  is  rendered 

X2-\-gi{xi)u  \ 

dt  \X2J  \Go{xi,X2) +  Gi{xi,X2)u-\-G2{xi)u^‘^W  ' 

where  Go{xi,X2)  =  go{xi,X2),  and  Gi  has  first  order  polynomial  dependence  on  X2.  The  Gi’s  are 
given  exjAicitly  by  equating  coefficients  in 

GoixuX2)  +  Gi{xi,X2)u  +  G2{xi)u^'^'^  =  go{xi,X2 -I  g\{xi)u) g2{xi-,X2  g\{x\)u)u 


+5r3(a:i)ut^^  - 


r  d 


dx\ 


{gi{xi)u) 


{X2  +  gi{xi)u). 


(12) 
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Proof:  It  is  easy  to  see  that  the  inverse  L-transformation  is  given  by 

/  a:i  \ 

\qj  \X2  +  gi{xi)uj 

Using  the  formulas  for  the  forward  and  inverse  L-transformations,  it  is  a  straightforward  calculation 
to  show  that  in  terms  of  (xi,  a:2)-coordinates  (3)  is  rendered  as  (11). 


Remark  3  In  the  case  that  our  system  is  given  by  (8),  the  L-transformation  is  the  ordinary  Leg¬ 
endre  transformation  provided  that  M{q)  =  I. 


It  is  not  difficult  to  see  that  (11)  inherits  a  polynomial  dependence  on  u  and  X2  which  is  similar  to 
(3).  The  main  feature  of  the  L-transformation  is  that  it  eliminates  the  dependence  on  the  derivative 
of  the  input  ii.  We  shall  next  see  that  if  ri(-)  is  a  high  frequency  periodic  input,  then  by  appropriately 
scaling  time  the  equations  (11)  may  be  rendered  in  a  form  to  which  averaging  theory  applies. 


Indeed,  suppose  u{-)  is  a  bounded  piecewise  continuous  R’^-valued  function  defined  on  (—00,00) 
which  is  periodic  of  fundamental  period  T  >  0.  Write 


Uj{t) 


E 


CjlcC  ^ 


k=  —  co 


and  let 


and 


The  corresponding  averaged  version  of  (11)  is 


(13) 


/  6+^i(6)w  y 

dt  V(f2/  \  Co(^l,  ^2)  +  Gi(^l,  1^2)^ -h  G2(<fi)cr^  / 

It  is  interesting  to  note  that  if  we  change  the  frequency  of  the  input  function  u(-),  the  form  of 
the  averaged  equations  remains  unchanged.  Specifically,  writing  u\(t)  =  u(\t)  we  see  that  u,\  has 
period  =  T/A  and  the  same  mean  u  =  u\dt  and  r.m.s.  value  u^^^dt.  In 

the  next  theorem,  we  shall  view  frequency  as  a  control  parameter. 


Suppose  (^1,^2) 


is  an  equilibrium  solution  of  (14).  Classical  averaging  theory  yields  the  following. 


Theorem  1  Let  u  (■)  an  -valued  bounded  piecewise  continuous  periodic  function  of  period 
T  >  0  whose  j-th  entry  is  given  by  (IS)  for  1  <  i  <  m.  Let  ((^i,<f2)  equilibrium  solution  of 

(14),  and  suppose  that  when  evaluated  at  (^1,(^2)  2n  x  2n-matrix 


(  iff" 

I  ^  ,  9^  2 

\  ^  96  “  ^  dii  ^ 


9Gn 

96 


I 

+ 
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has  all  its  eigenvalues  in  the  left  half-plane.  Then  for  every  e  >  0  there  is  a  6  >  0  and  Aq  such  that 
for  all  A  >  Ao,  if  (11)  is  forced  by  the  input  ux{t)  =  u{Xt),  we  find  the  solution  {xi{t),X2{t))  of  (11) 
is  such  that  ||a;(f)  -  <  e  provided  ||a;(0)  -  <  <5. 


Proof:  Given  e,  6,  and  A  as  in  the  hypothesis  of  the  theorem,  we  shall  write  both  (11)  and  (14)  in 
the  time-scale  r  =  t/Tx  where  Tx  =  j  is  the  fundamental  period  of  ua(-)-  Letting  2/(r)  =  x{t),  (11) 
becomes 


dr  \y2 

(Cf.  (9).)  Similarly,  letting  ((t)  =  ((t),  (14)  becomes 


\  _rp  f  y2  +  giiyi)ux{TTx) 

1/  ^  \Go{yi,y2) Gi{yi,y2)ux{T'Tx)  +  G2{yi)ux{TTx)^^^ 


±  (Cl 

dt  \C2 


=  Tx 


C2  +  </i(Ci)^ 

Co(Cl)  C2)  +  Gi(Cl,  C2)ii  +  G2{Ci)^^ 


(15) 


(16) 


Note  that  {i(r)  =  ux{tTx)  is  a  periodic  function  of  period  27r  independent  of  A. 


Now  (15)  is  a  nonautonomous  equation  to  which  classical  averaging  theory  may  be  applied.  (Specif¬ 
ically,  see  [5],  p.  497.)  Under  our  hypothesis,  the  eigenvalues  of  the  Jacobian  of  the  right  hand 
side  vector  field  (16),  evaluated  at  (^J,^”),  all  have  negative  real  parts.  Hence  for  all  A  suffi¬ 
ciently  large  (and  Tx  =  T/\  correspondingly  small),  there  is  a  unique  periodic  (period  27r)  solution 
(?/i(t),  t/2('r))  to  (15)  defined  on  the  entire  infinite  interval  (—00,00)  such  that  ||j/(r)  —  i^°l|oo  <  vW 
for  -00  <  T  <  00  with  7]{\)  0  as  A  ->  00.  (See  [5].)  Moreover,  if  y{T)  is  any  solution  to 

(15)  which  is  sufficiently  close  to  C°  at  some  time  r  =  tq,  then  ||y(r)  —  ?/(r)||  — >  0  exponentially  as 
T  00.  There  is  no  loss  of  generality  in  assuming  tq  =  0. 


To  complete  the  proof  of  the  theorem,  let  A  be  sufficiently  large  that  77(A)  <  Moreover,  let 
Ti  be  such  that  for  all  t  >  tj  we  have  ||y(r)  —  y(r)||  <  e/2.  Let  5  >  0  be  so  small  that  if 
||?/(0)  -  y{0)\\  <  5,  then  ||7/(t)  -  y(T)ll  <  e/2  on  0  <  r  <  Tj.  Now  A  has  been  chosen  large  enough 
to  ensure  that  77(A)  <  f,  and  hence  |j7/(r)  -  C'^Woo  <  f  •  Thus  if  ||y(0)  -  ^°||  <  f ,  we  have 

ll2/(0)-y(0)ll  <  ||7/(0)-^°||  +  ||y(0)-^°ll 
S  6 
^  2^2 

=  6.  (IT) 

Hence  ||y(T)  -  7/(t)||  <  f  for  all  t  >  0,  and  thus  lly(T)  -  <  e.  Since  x  and  y  are  related  by  a 

change  of  time  scale,  this  proves  the  theorem.  ^ 


Remark  4  Under  the  hypothesis  of  the  theorem,  we  have  proved  that  solutions  of  (11)  will  execute 
motions  confined  to  a  neighborhood  of  (^?,  ^2)-  Moreover,  as  the  frequency  of  the  input  function  u(-) 
becomes  larger,  this  neighborhood  to  which  the  motions  are  confined  shrinks.  We  may  relate  these 
observations  to  motions  of  (3).  The  variables  {q,q)  are  related  to  (xi,a;2)  by  means  of  the  inverse 
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L-transformation,  and  it  is  not  difficult  to  see  that  trajectories  expressed  in  these  coordinates  may 
be  expected  to  execute  motions  in  a  neighborhood  of  the  point 

Since  the  trajectory  {q,q)  is  related  to  the  trajectory  {xi,X2)  by  means  of  this  time- varying  co¬ 
ordinate  change,  the  size  of  the  neighborhood  to  which  the  (g,  ^j-motions  are  confined  does  not 
necessarily  shrink  as  the  frequency  of  the  forcing  is  increased. 


Remark  5  It  is  important  to  note  that  this  theorem  does  not  actually  deal  with  stability  of  equi¬ 
libria  for  systems  of  the  form  (11)  since  (<^1,(2)  needs  only  to  be  an  equilibrium  (stable)  for  the 
averaged  system  (14).  It  is  easy  to  find  examples  (See  below!)  where  {xi,X2)  =  (i^i,<^2)  not  an 
equilibrium  of  (11)  and  yet  (11),  under  the  hypotheses  of  the  theorem,  executes  motions  confined 
to  a  neighborhood  of  (i^?,^®)- 


Example  5  (Example  1  reprise) 


The  averaged  equations  (14)  in  this  case  take  the  form 


d 

dt 


Clearly  the  hypotheses  of  the  theorem  are  not  met,  and  no  conclusion  can  be  drawn.  A  possible 
way  to  further  analyze  this  system  is  to  pass  to  higher  order  averaging  approximations  as  discussed 
in  [5].  Equivalently,  we  may  also  pass  to  the  representation  of  the  system  given  in  Example  2. 


Example  6  (Example  2  reprise)  For  this  representation,  the  averaged  equations  (14)  take  the 
form 

/C  \  6  +  hsinC  A 

df  \  6  /  \—u  cos  C  •  <^2  —  <7^  sin  C  cos  C  J  ' 

Taking  x{t)  =  cos  at,  we  have  u{t)  =  sin  of,  and  ii  =  0,  =  |.  (>^1,1^2)  =  (O5O)  is  clearly  an 
equilibrium  solution  of  the  averaged  equations.  The  linearized  version  of  this  differential  equation 
at  this  equilibrium  has  eigenvalues  on  the  imaginary  axis,  and  hence  the  hypothesis  of  Theorem  1 
is  again  not  satisfied.  We  may  proceed,  however,  by  appealing  to  physical  considerations. 


Returning  to  the  description  of  the  cart-pendulum  system  in  the  introduction,  we  note  that  if  a 
rate-dependent  damping  term,  dO,  is  included  in  the  model,  (2)  becomes 

iO  +  dO  —  X  smO  =  0. 

In  terms  of  the  general  representation  (3),  ii  =  x,  go{0,d)  =  —d6,  gi{0)  =  sin^,  g2  =  0,  and  g^  =  0. 
The  corresponding  averaged  equations  are 

6  +  «sinC  \ 

V6/  \-d{(2  +  usm(i)  -ucos^i  •  ^  sin ,^1  cos ^1  /  ' 


dt 
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Linearizing  about  the  equilibrium  =  (0^0)  yields  the  variational  equation 

u  1 

dt\8^2J  K-du-a"^  -d-u)  \S^2)  ^ 

the  coefficient  matrix  of  which  has  characteristic  polynomial  p{s)  =  s'^  +  ds a'^  —  .  Noting  that 

(7^  —  >  0,  the  eigenvalues  are  seen  to  lie  in  the  open  left  half  plane.  Hence,  by  our  theorem,  the 

forced  system  can  be  expected  to  execute  stable  motions  in  a  neighborhood  of  the  origin. 


This  observation  is  consistent  with  what  we  would  expect  to  observe  if  viscous  damping  were 
added  to  the  dynamics  of  a  Lyapunov  stable  system.  It  remains  to  show,  however,  that  the  origin 
[0^0)  -  (0,0)  is  Lyapunov  stable  for  the  periodically  forced  system  (2).  While  the  above  analysis 
does  not  directly  imply  this,  it  suggests  a  general  approach  to  stability  analysis  for  conservative 
systems.  The  idea  is  to  add  Rayleigh  dissipation  (along  the  lines  suggested,  say,  in  [1])  and  then 
assess  the  spectral  stability  of  the  resulting  averaged  system.  Similar  ideas  have  recently  been 
explored  in  [4]  in  connection  with  the  energy-momentum  method.  We  shall  provide  a  brief  general 
sketch. 

For  systems  of  the  form  (3)  dissipation  enters  through  the  terms  gl{q)q  (see  (4)).  Indeed,  we  shall 
assume  that  gl{q)  is  a  symmetric  negative  definite  matrix  so  that  —q^gl{q)q  is  the  usual  Rayleigh 
dissipation  function.  We  shall  be  interested  in  comparing  systems  of  the  form  (3)  which  differ  from 
one  another  precisely  according  to  whether  the  Rayleigh  dissipation  function  is  absent  or  present. 
Accordingly,  we  rewrite  the  equation  of  interest,  (3),  in  the  form 

q  =  9l{q)  +  ^9l{q)q  +  9l{q)q^^^  +  9i{q)'d  +  92{q,  +  93,{qW^^^  (is) 

where  e  =  0, 1  depending  on  whether  damping  is  absent  or  present.  Applying  the  L-transformation 
(10)  and  averaging,  we  obtain  from  (18) 

—  6  +  S'i(6)ll  QQ'I 

dt  V6/  \^9o{^i){^2  +  9iUi)d)  +  Go(6i'^2)  +  ’ 

where  the  G^s  are  related  to  the  gds  via  equation  (12)  with  e  =  0  in  go{^i).  Suppose 
an  equilibrium  for  this  differential  equation;  i.e.  ^2  +  9i{^i)d  —  0  ^^i^d  Go(^i5^2)  +  + 

G2(^?)(T^  =  0.  Then  we  have  the  following. 

Theorem  2  For  the  conservative  version  of  (19)  (i.e.  when  e  =  0),  let  us  suppose  that  for  some 
admissible  control  input  u{-)  (i.e.  for  some  bounded  piecewise  analytic  u)  and  corresponding  u,  a 
there  is  a  function  E  of  the  form 

Eifuf2;u,a)  =  ^116  +  5'i(6)w||^  + 1^(6; 

which  (i)  has  a  strict  local  minimum  at  the  equilibrium  (if?, ^2)/  (ii)  invariant  under  the 
motion  of  (19).  Then  for  the  dissipative  version  of  (19)  (i.e.  when  we  set  e  =  1)  the  equilibrium 
(f?,f?)  is  asymptotically  stable. 
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Proof:  Differentiating  E  along  trajectories  of  (19)  we  obtain 

E  =  «(6  +9.«l)iif9j({i)(f2  +  Sl(6)6). 

Here  E  serves  as  a  Lyapunov  function,  and  the  right  hand  side  of  the  equation  is  nonpositive  and 
has  a  strict  local  maximum  value  of  0  at  (^1,^2)-  By  Lyapunov’s  second  method,  this  equilibrium 
solution  is  asymptotically  stable.  □ 


Remark  6  The  function  E  plays  the  role  of  total  energy  for  the  averaged  system  (19).  This  role 
may  be  emphasized  by  applying  the  averaged  inverse  L-transformation 

it)  ^  {^2+tui)u) 

and  expressing  E  in  terms  of  the  transformed  variables  {q,q): 

E(q,^-,u,a)  =  +  V{q-,ii,a). 


Corollary  1  For  the  conservative  version  of  (19)  (i.e.  when  e  —  0),  suppose  that  for  some  admissi¬ 
ble  control  input  u(-)  and  corresponding  u,  a  there  is  a  function  E{^i  cr)  as  in  Theorem  2  which 

has  a  strict  local  minimum  at  the  equilibrium  (^?,^2)  which  is  invariant  under  (19).  Then  for 
every  q  >  0  there  is  a  6  >  0  and  Aq  such  that  for  all  X  >  Xq  if  the  corresponding  dissipative  version 
of  (11)  (derived  by  applying  the  L-transformation  to  (18)  with  e  =  1)  is  forced  by  u\{t)  =  u{Xt),  the 
solution  x{t)  =  {xi{t),X2{t))  is  such  that  ||a:(t)  —  <  9  provided  ||a:(0)  —  (f°||  <  8. 


Proof:  This  corollary  follows  directly  from  Theorems  1  and  2  and  standard  Lyapunov  stability 
theory.  □ 


Remark  7  In  Lagrangian  and  Hamiltonian  mechanics,  the  most  obvious  choice  of  function  which  is 
invariant  under  the  dynamics  of  the  system  is  the  total  energy.  Thus  to  analyze  mechanical  systems 
with  inputs  of  the  form  (3),  it  is  of  interest  to  ask  whether  there  is  any  systematic  procedure  for 
relating  the  corresponding  averaged  system  (14)  to  a  Hamiltonian.  (Because  inputs  are  present,  the 
usual  total  energy  is  not  conserved.)  In  the  next  section,  it  will  be  shown  that  the  reduced-order 
symmetric  Lagrangian  control  systems  of  Example  4  retain  a  structure  which  may  be  related  via 
the  Legendre  transformation  to  a  Hamiltonian  system.  Averaging  this  system  leads  to  another 
Hamiltonian  system,  so  that  stability  analysis  along  the  lines  we  have  indicated  may  be  carried 
out.  Moreover,  the  analysis  will  be  considerably  simplified  for  this  class  of  systems  because,  as  we 
shall  show,  there  is  a  natural  block-diagonalizing  coordinate  transformation  such  that  in  the  new 
coordinates  the  total  averaged  energy  (=the  averaged  Hamiltonian)  may  be  expressed  as  the  sum  of 
a  positive  definite  averaged  kinetic  energy  function  and  an  averaged  potential  function.  Equilibrium 
points  of  the  averaged  system  will  be  found  as  critical  points  of  the  averaged  potential,  and  the 
stable  equilibria  of  interest  in  our  analysis  will  be  local  minima  of  the  averaged  potential. 
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4.  Energy  Methods  for  Stability  and  the  Averaged 

Potential 


Consider  a  symmetric  Lagrangian  control  system  of  the  form  (6)  where  qi  is  a  cyclic  variable  as 
defined  in  Example  4.  Write  v  -  qi.  As  in  the  preceding  section,  we  view  u  as  a  control,  and  note 
that  (8)  arises  from  applying  the  Euler-Lagrange  operator  to  the  reduced  Lagrangian 

L{q,  q-,  u)  =  ^fM{q)q  +  fA{q)v  -  Va{q\  u),  (20) 

where  Va{q\v)  =  V{q)  —  \v'^m{q)v  is  the  augmented  potential. 

To  render  the  equations  of  motion  in  a  form  suitable  for  averaging,  we  apply  the  Legendre  trans¬ 
formation  Lf  (5,  p;  v)  =  p  ■  q  —  L{q,  q\  v)  to  obtain 

H(q,p-,v)  -  ^{p  -  vAfM~'^{p  -  vA)  +  Va{v,q).  (21) 

The  equations  of  motion  are  written  in  the  usual  way  as 

q  =  M~^{p  —  vA)  (22) 

p  =  -^\^-{p  -  vAfM-\p  -  vA)  -b  K(u, q)].  (23) 

oq  I 

To  pursue  averaging  theory,  we  replace  all  coefficients  in  (22)-(23)  by  their  time-averages.  Assuming 
?;(•)  is  bounded,  piecewise  continuous,  and  periodic  of  period  T  >  0,  we  may  write 

00 

^^(0=  E  (24) 

Ar=  — CO 

(Cf.  (13).)  (22)-(23)  contain  terms  of  order  not  greater  than  two  in  u,  and  averaging  the  coefficients 
we  obtain 


Proposition  2  Suppose  is  given  by  (24)-  Then  if  all  coefficients  in  (22)-(23)  are  replaced  by 
their  time  averages,  the  resulting  averaged  system  is  Hamiltonian  with  corresponding  Hamiltonian 
function 

H{q,p)  =  ^(p  -  A{q)vf  M{q)~^{p  -  A(?)u)  +  VA{q), 

where 

VA{q)  =  V{q)  +  ^(S(g)  -  v'^A{qfM{q)-^A{q)v) 
is  the  averaged  potential  with 

1  fT 

”  =  fi„ 


S(?)  =  t  f  v(tf(A(qfM(q)  'A(q)-m(q))v{t)dt. 
1  Jo 


and 
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□ 


The  averaged  Hamiltonian  H{q,p)  may  be  interpreted  as  an  averaged  total  energy.  Writing  this  is 
terms  of  the  variables  q  and  q,  where  q  =  M~^{p  —  Av)  the  averaged  energy  is 

E{q,q)  =  ^fM{q)q+VA{q). 

Equilibrium  solutions  to  the  averaged  Hamiltonian  system 

.  dH  .  dli 

will  be  of  the  form  (9,  q)  =  ((/o?  0).  We  have  the  following: 


Proposition  3  Let  (q^p)  —  (qo^Po)  be  an  equilihrhuii  solution  of  the  averaged  Hamiltonian  system 
(25).  {qo>)Po)  is  a  strict  local  minimum  of  the  total  averaged  energy  H{q.^p)  if  and  only  if  qo  is  a 
strict  local  minimum  of  the  averaged  potential  VA{q). 


Proof:  The  corresponding  equilibrium  in  (g,  ^)-coordinates  is  (50,0).  (^o^Po)  is  a  strict  local  min¬ 
imum -of  H{q.^p)  if  and  only  if  (go?0)  is  a  strict  local  minimum  of  E{q.,q).  But  the  Hessian  of  E 
evaluated  at  (^0?  0)  has  block  diagonal  form 


( 0  \ 

\  0  M(,„)p 

The  proposition  follows  from  this  observation. 


□ 


To  discuss  stability,  we  introduce  the  Rayleigh  dissipation  function  q^D{q)q  where  for  each  D{q) 
is  an  n  X  n  positive  definite  matrix.  Starting  from  the  reduced  Lagrangian  (20),  we  are  interested 
in  the  behavior  of  the  system 


d  dL 
dt  dq 


(26) 


The  Hamiltonian  form  of  (26)  is  given  by 

.  dH 

B  H 

p= -  D{q)M{q)~^{p- A{q)v),  (27) 

where  H  is  as  in  (21).  Averaging  these  equations,  we  obtain  the  system 
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p=-—-D{q)Miq)-\p-A{q)v),  (28) 

dq 

where  the  averaged  Hamiltonian  H  is  given  in  Proposition  2.  The  role  of  the  averaged  potential  in 
assessing  asymptotic  stability  is  summarized  as  follows. 


Theorem  3  Consider  a  Lagrangian  control  system  (20)  with  dissipation  entering  through  a  quadratic 
dissipation  function  q^D{q)q  as  described  above.  Let  v{-)  be  given  by  (24),  let  vx{t)  =  v{\t),  and  for 
each  X,  consider  the  averaged  potential  VA{q;X)  corresponding  to  vx{-)  as  in  Proposition  2.  Suppose 
further  that  there  is  a  Aq  and  a  one-parameter  family  of  configurations  qx  depending  smoothly  on  A 
with  the  property  that  for  each  A  >  Aq  the  corresponding  qx  is  a  strict  local  minimum  o/VA((jf;A). 
Then  there  is  a  Ai  >  Aq  such  that  for  each  e  >  0  and  A  >  Ai  there  is  a  8  =  6{e,  A)  with  the  property 
that  if  (27)  is  forced  by  vx{-),  the  solution  {q{t),q{t))  satisfies 

||(<?,9)  -  (<ZA,0)||oo  <  e 


provided 


Proof:  Let  e  >  0  be  given,  and  let  u(-)  be  defined  by  (24)  with  corresponding  ua(-)  as  in  the 
statement  of  the  theorem.  Let  Hx  be  the  corresponding  averaged  Hamiltonian.  If  we  let  g,  q  = 
M{q)~^{p  —  A{q)v)  be  the  corresponding  Lagrangian  coordinates  of  the  averaged  system,  Hx  may 
be  rewritten  in  these  coordinates  as 

E\{q,q)  =  ]^q^M{q)q  +  VA{q\X). 


In  ((?,  g)-coordinates,  (28)  is  rendered  as 

f{M(q)q)  -  ^j\fM(q)q  +  VAq\  A)1  +  D(q)q  =  0.  (29) 


Differentiating  Ex{q,q)  along  trajectories  of  this  system,  we  find 


dEx 

dt 


-q'^D{q)q. 


D{q)  is  positive  definite  for  each  q,  and  since  {q,q)  =  (^^,0)  is  a  strict  local  minimum  of  Ex,  it 
follows  that  {qx,0)  is  an  asymptotically  stable  equilibrium  of  the  averaged  system  (29). 


It  is  now  more  or  less  straightforward  to  apply  Theorem  1  to  determine  the  behavior  of  (26).  To 
do  this,  apply  the  L-transformation  to  the  variables  {q,q).  Now  the  L-transformation  is  related  to 
the  Legendre  transformation  by  a  left  multiplication  by  d\&g{I,M{q)~^).  Hence  the  asymptotically 
stable  equilibrium  solution,  (^a^O),  of  (29)  corresponds  to  an  asymptotically  stable  equilibrium 
solution  of  (14).  It  thus  follows  that  (26)  executes  motions  in  a  neighborhood  of  {q\,0).  □ 
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Figure  2:  A  cart-pendulum  system  on  an  inclined  track. 


Remark  8  The  form  of  the  averaged  system  (29)  may  be  determined  by  examining  the  form  of  the 
system  (26)  from  which  it  is  derived.  Specifically,  in  passing  from  (26)  to  (29),  terms  of  the  form 


E 


+  vj{ 


i=i 


dq 


(.30) 


are  replaced  by 

^  ~lm  -  -  mW)ii|. 


(31) 


Example  7  (Example  6  reprise.)  It  is  of  interest  to  consider  an  extension  of  the  cart-pendulum 
problem  which  serves  to  illustrate  the  power  of  the  theory  developed  so  far.  Consider  a  cart  to  which 
there  is  a  pendulum  attached  moving  along  a  track  in  a  uniform  gravitational  field.  Suppose  the 
track  is  inclined  at  an  angle  a  with  respect  to  the  horizontal.  As  in  the  introduction,  the  pendulum 
consists  of  a  point-mass  attached  by  a  massless  link  of  length  t  to  the  cart  which  has  mass  iric.  The 
coupled  equations  of  motion  for  this  system  are 

(rUc  -f  mp)r  -b  rUpi cos{6  —  a)6  —  rripl sin(0  —  a)9'^  -f  {rric  -b  mp)g  sin  a  —  u,  (32) 

and 

-b  r  cos{9  —  Of)  -b  5'  sin  ^  =  0,  (33) 

which  reflect  the  assumption  that  control  of  this  system  is  by  pushing  the  cart  along  its  inclined 
track.  As  in  our  earlier  analysis,  we  may  assume  that  the  motion  r{t)  can  be  specified,  and  taking 
v{t)  =  r{t),  we  find  that  (33)  is  derivable  from  the  reduced  Lagrangian 

L{9, 9;  v)  =  -i9^  +  v9  cos(9  —  a)  +  g  cos  9.  (34) 

Assume  v{t)  is  given  as  in  (24).  Referring  to  the  formulas  in  Proposition  2,  we  write  the  averaged 
potential  corresponding  to  L: 


24 


BAILLIEUL 


where  ^ 

f  Jo 

and  y 

Stable  equilibria  of  the  averaged  system  correspond  to  local  minima  of  Va-  One  interesting  special 
case  occurs  when  g  =  0  (i.e.  the  gravitational  field  is  removed),  o;  =  0,  and  v{t)  ~  sinut.  This 
corresponds  to  the  example  treated  in  the  introduction.  Then  u  =  0  and  cr^  =  5,  and  0  =  ±|  are 
local  minima  of  Va-  If  there  is  any  dissipation  in  the  hinge  on  the  pendulum,  then  according  to 
Theorem  3,  the  pendulum  may  be  expected  to  execute  motions  confined  to  a  neighborhood  of  |  or 
Another  interesting  special  case  is  to  take  a  =  |  with  nonzero  g  and  r{t)  =  -cos(wt).  (The 
cart  executes  a  simple  periodic  vertical  motion  along  its  track.)  Then  v{t)  =  h  =  0,  and 

cr^  =  ^  so  that  the  critical  points  of  the  averaged  potential  are  solutions  of 

g  sin  ^  ^  0  =  0. 

This  equation  has  solutions  0  =  0,  tt  for  all  values  of  the  parameters,  and  when  >  2gi  there  are 
also  solutions  0  =  arccos(-2^).  shall  restrict  0  G  (-7r,;r].)  Evaluating  the  second  derivative 
of  Va  at  each  of  these  critical  points,  we  find  that  0  =  0  is  always  a  local  minimum,  and  0  =  tt  is  a 
local  maximum  when  <  2gi,  but  it  becomes  a  local  minimum  when  >  2gi.  It  follows  from 
Theorem  3,  that  there  is  an  uq  >  2g£  such  that  for  all  u>  >  uq  the  pendulum  can  be  made  to  execute 
motions  in  a  neighborhood  of  the  vertical  configuration.  (Cf.  [5],  p.  408.)  A  similar  analysis  can 
be  carried  out  for  other  values  of  the  parameter  a.  We  refer  the  reader  to  [15]  for  further  details. 


5.  The  Geometry  of  Averaged  Mechanical  Systems 


If  we  assume  the  input  function  u(-)  discussed  in  the  preceding  section  has  the  form  (24),  then 
we  see  that  the  averaged  potential  Va  is  invariant  under  a  scaling  of  frequency  v{t)  — >  v{Xt)  for 
A  0.  Hence  the  form  of  Va  contains  only  partial  information  regarding  the  stability  of  the  system, 
since  the  hypothesis  of  Theorem  3  requires  both  that  A  is  large  enough  so  that  the  critical  point 
of  VA  is  a  strict  local  minimum,  and  that  A  is  large  enough  so  that  classical  averaging  theory 
applies.  Nevertheless,  further  analysis  may  proceed  based  on  knowledge  of  Va  alone.  Based  on  the 
conclusions  of  Theorem  3,  we  might  conjecture  that  if  the  frequency  of  the  forcing  is  sufficiently 
high  and  if  qo  is  a  strict  local  minimum  of  Va,  then 


dt  dq  dq 


(35) 


will  execute  motions  in  a  neighborhood  of  {q,  q)  =  {qo,  0).  An  essential  assumption  in  proving  The¬ 
orem  3  was  that  dissipation  was  present,  and  we  are  now  interested  in  investigating  circumstances 
under  which  the  result  also  holds  for  the  conservative  system.  The  required  analysis  must  be  car¬ 
ried  out  along  different  lines  since  there  is  no  way  to  take  advantage  of  asymptotic  stability  of  a 
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linearization.  We  shall  show,  for  a  certain  prototypical  system,  that  in  the  absence  of  dissipation, 
stability  will  be  the  effect  of  a  nonholonomic  response  of  the  system  to  the  periodic  forcing.  To 
pursue  our  analysis,  we  linearize  (35)  about  the  point  (^OiO):  Substituting  q  =  qQ-\-  8q  into  (35)  and 
retaining  terms  up  to  first  order  in  Sq,  we  obtain 

M{qo)Sq+^{qQ,v)  +  ^^{qo,v)8q 

+  +  ^(9o)^9)  +  ^  )4  =  0.  (36) 

Two  features  of  this  system  are  of  interest:  (i)  It  is  jointly  bilinear  in  the  inputs  u,  v  and  states  q, 
4,  and  (ii)  the  linearization  is  about  a  point  which  may  or  may  not  be  an  equilibrium  of  the  system. 
When  the  point  qo  is  not  an  equilibrium  of  the  forced  system,  V.  Solo  ([16])  has  suggested  the  term 
“hovering  motion”  replace  the  term  “stable  motion.”  We  shall  pursue  the  analysis  of  both  types  of 
motion  under  the  general  heading  of  “stability.”  The  advantage  of  studying  the  system  (36)  is  that 
stability  analysis  may  be  pursued  using  powerful  linear  tools  such  as  Floquet  theory.  This  may  be 
illustrated  in  terms  of  the  cart-pendulum  problem. 


Example  8  (Example  7  reprise)  Recall  (34) 

L{9, 6]  v)  —  -iO^  +  v6  cos[0  —  a)  +  g  cos  6. 


Specializing  the  terms  in  (36),  we  have  M{0)  =  f,  A{6)  —  cos{6  —  a),  and  Va{9)  =  —geos  0.  Hence, 
we  may  write  the  variational  equation  corresponding  to  (33)  with  v  —  r  as 


88e\ 

dt  \80  ) 


(  “  ■ 

V  -i  ^g  cos  9o+  j  sin(^o  -  ct) 


0 

jcos{9o  -  a) 


(37) 


Pursuing  an  analysis  of  the  second  case  considered  in  Example  6,  we  set  o  =  ^  and  consider  the 
critical  point  9o  =  of  the  averaged  potential.  The  above  equations  specialize  to 


d_  f89\  _  f  0  1\  f89\ 

dt\89)-[s^  0)\89)- 

In  general,  if  we  have  a  linear  differential  equation  with  periodic  coefficients  of  the  form 


(38) 


x{t)  —  F{t)x{t),  (39) 

where  the  coefficient  matrix  has  entries  which  are  bounded,  piecewise  analytic,  and  periodic  of 
period  T  >  0,  then  Floquet  theory  implies  that  the  trivial  solution  to  (39)  will  be  (Lyapunov) 
stable  if  all  eigenvalues  of  the  fundamental  matrix  $(7",  0)  are  distinct  and  lie  on  the  unit  circle. 
When  this  occurs,  we  have  the  following 


Proposition  4  Let  $(t,0)  denote  the  fundamental  matrix  solution  to  (29),  and  suppose  $(T’,  0) 
(T  >  0  is  the  fundamental  period  of  v)  has  distinct  eigenvalues  lying  on  the  unit  circle.  The 
solution  to  (39)  sampled  every  T  units  of  time  evolves  on  an  ellipsoid  x^Mx  =  1. 
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Proof:  If  all  eigenvalues  of  $(r,0)  are  distinct  and  lie  on  the  unit  circle,  there  is  a  nonsingular 
matrix  P  such  that  A  =  P“^$(r,0)P  is  block  diagonal  with  O’s  and  2  x  2  blocks  of  the  form 

f  cos  -  sin  ^  diagonal.  Now  AA^  =  I,  and  from  this  it  follows  that  $(r,  0)^M$(r,  0)  = 

V  sin  <p  cos  <p  J 

M  where  M~^  =  PP^.  With  the  appropriate  scaling  of  M,  the  proposition  follows.  □ 

Suppose  we  introduce  the  phase-shift  v(t)  i— >  v{t  -f  e).  Then  it  is  easy  to  show  that  the  time-T 
invariant  ellipsoid  is  defined  by  the  quadratic  form  =  $(0,  e)^M$(0,  c).  The  union  of  the  family 
of  invariant  ellipsoids,  \j{x'^MeX  =  1  :  0  <  e  <  T},  provides  an  approximation  of  the  region  in 
the  state-space  through  which  the  trajectories  of  the  forced  system  (39)  pass.  Indeed,  if  a;(0)  =  xq 
and  x^Mxo  ==  1,  then  since  x{s)'^Msx{s)  =  x^Mxq,  the  solution  to  (39)  with  initial  conditions 
a;(0)  =  Xq  lies  in  this  union  of  ellipsoids. 


A  general  theory  of  how  well  this  set  describes  the  set  of  points  through  which  the  solution  passes 
does  not  presently  exist.  The  cart-pendulum  problem  illustrates  what  is  involved,  however. 


Example  9  (Example  8  reprise)  Consider  the  linearized  dynamics  for  the  cart-pendulum  problem 
(38).  Letting  u{t)  =  -{g  -f  v{t))/£  and  y  =  (?/i,y2)  =  {SO,  60),  this  represents  the  well-studied 
Hill’s  equation.  (See  e.g.  [11].)  We  shall  conclude  this  section  with  a  discussion  of  the  relationship 
between  the  geometry  of  solutions  of  (37)  and  the  averaged  potential  introduced  in  the  preceding 
section.  Suppose  the  input  u  is  a  square  wave 


u(t) 


-7 

[  u{t  —  2h) 


ii  0  <  t  <  h] 
if  h  <  t  <  2/i; 
if  2h  <  t. 


{j3,  7  >  0).  The  fundamental  matrix  $(t,0)  has  determinant  1  for  all  t,  and  at  time  T  =  2h  may 
be  explicitly  written  as 


^2h,0) 


/  cosy/^h  -^s'my/^h\  /  coshy^/i  ^sinh^/7/i\ 
[-y/^Smy/^h  COSy/^h  )  [^siuh  COsh  )' 


As  our  previous  analysis  suggests,  the  solution  (j/i,j/2)  =  (0)0)  will  be  stable  precisely  when  the 
eigenvalues  of  $(2/j,0)  are  on  the  unit  circle.  The  characteristic  polynomial  of  $(2/i,0)  is 


p{s)  =  5^  -  /i(/?,7,/i)s  +  1 


where  p{(3,  'f,h)  =  2  cos  y/^h  cosh  ^h  +  {^ -  y/^)  sin  y/^h  sinh  y/^h.  The  roots  of  this  polynomial 
will  be  on  the  unit  circle  precisely  when  |/i(/?,7, /i)|  <  2.  When  /?  =  7  =  1,  for  example,  one  sees 
that  there  is  a  sequence  of  roots  of  the  equation  )u(l,  1,  h)  =  2:  /lo  =  0,hi  =  1.8751,  /i2  =  4.73004, . . . 
with  — >  00  as  n  ^  00.  In  the  interval  (/i2A:)  ^2fc+i))  the  trivial  solution  to  Hill’s  equation  ((38) 

will  be  stable,  while  for  h  €  (/i2fc+i)  h2k+‘2),  the  solution  will  be  unstable. 


Applying  the  analysis  of  the  preceding  section  to  this  system.  Theorem  3  implies  that  for  all  suffi¬ 
ciently  small  h  >  0,  the  origin  will  be  a  stable  equilibrium  for  (38).  A  more  refined  picture  emerges 
from  an  analysis  of  the  invariant  ellipsoids. 
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Proposition  5  Given  a  2  x  2  matrix  $  =  with  the  property  that  det$  =  1,  there  is  a 

symmetric  2x2  solution  M  =  {  )  to  the  matrix  equation 

\mi2  m22/ 

=  M 


given  explicitly  by 

(  ”^11  \ 
mu  = 
\  ”*22 

This  solution  is  unique  to  within  a  scale  factor. 


(  —2(j)2l  \ 
fll  —  <^22 

\  ‘2.fu  ) 


Proof:  A  direct  computation,  using  the  hypothesis  det  $  =  1,  verifies  that  the  solution  we  have 
proposed  does  indeed  satisfy  the  equation.  To  prove  that  this  solution  is  unique  up  to  a  scale  factor, 
it  is  only  necessary  to  show  that  the  linear  mapping  M  i->  has  rank  2.  To  do  this,  we  shall 

prove  that  the  images  are  linearly  independent.  Suppose  we  have  scalars 

a  and  h  such  that 

/  <t>\\  4>n<i>\2\,,(  <?i>21  <r’21<?!>22\  _  /O  0\ 

Kfllfu  fu  /  \f2\f22  ^22  /  VO  0/ 

It  is  straightforward  to  show  that  this  equation  is  equivalent  to 

/  <t>\\  4>2\  \  f  afii  A  /  ?i’ii  <t>2\  V  f  afu  \  /  0  V 

\4>12  4>22J\bf2\)  \fl2  <l>22  J  \  b(f>22  )  VO/' 

Since  det  $  =  1,  this  is  equivalent  to 

/  a4>u  \  ^  ( a(j)u  \  _  /OV 
V  b(j)2\  )  V  h(j}22  J  VO/ 

Since  we  cannot  have  (/in  =  <^21  =  0,  we  must  have  a  =  0.  Similarly,  6  =  0,  proving  the  rank  of  the 
linear  mapping  is  2,  as  claimed.  This  proves  the  proposition.  □ 


Remark  9  Given  $(2/i,0)  as  in  Example  9,  let  M  be  the  corresponding  invariant  quadratic  form 
specified  by  Proposition  5.  The  condition  that  M  may,  after  a  multiplication  by  a  suitable  scalar, 
be  taken  to  be  positive  definite  is  that  m\im22  —  >  0.  Note  that  this  is  equivalent  to  the 

condition  that  $(26,0)  has  distinct  eigenvalues  on  the  unit  circle:  i.e.  to  {<f>ii  +  </>22)^  —  4  <  0.  □ 


Given  a  square  wave  input  u(-)  as  prescribed  above,  we  write  the  corresponding  solution  to 

$(26,0)^M$(26,0)  =  M, 

and  normalize  M  with  respect  to  the  initial  conditions  such  that  y{0)^ My{0)  =  1.  If  M  is  positive 
definite,  solutions  (7/1 , 7/2)  =  {bO,  69)  of  (38)  which  satisfy  this  equation  evolve  on  the  union  of  ellipses 
G{y^M(^y  =  1  :  0  <  e  <  T}.  A  typical  motion  is  shown  in  Figure  3,  where  the  ellipses  y^Moy  =  1 
and  y"^ Mhy  =  1  are  also  depicted. 
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Figure  3:  Phase  portrait  of  solution  to  Hill’s  equation  with  two  members  of  the  family  of  invariant 
ellipsoids.  Forcing  is  by  square  wave  function  u{t)  as  described  in  the  text  with  parameters  a  = 
/?  =  1  and  h  =  0.3. 


Our  final  objective  will  be  to  return  to  the  question  of  how  our  results  vary  as  the  the  frequency 
of  the  forcing  is  changed.  This  will  provide  an  important  complement  to  the  discussion  in  the 
previous  section  where  it  was  noted  that  the  averaged  potential  Va  is  invariant  under  a  change  in 
the  frequency  of  the  forcing.  Returning  to  (38),  it  will  be  useful  to  now  assume  v{t)  (instead  of 
u{t))  is  the  square  wave 

(  —7  if  0  <  t  <  /j; 

v{t)  =  \  13  \{h<t<  2h-,  (40) 

[  v{t  —  2h)  if  2/1  <  t, 

(y3,  7  >  0).  To  study  the  effect  of  scaling  the  frequency,  we  define  v\{t)  =  v{\t)  so  that 


v\it) 


'  -A7  if  0  <  t 

<  if  f  t  <  ^;  . 

^vx{t-2hlX)  ift>f. 


From  Proposition  2,  we  write  the  averaged  potential  for  the  cart-pendulum  system  (34)  in  the  case 
0  =  1  and  the  input  has  7  =  /?: 

Va^O)  =  -S'  cos  6>  +  ^  sin^  0  •  (o^  -  u^). 


A  straightforward  computation  to  evaluate  u,  a  yields 

Va{&)  = -g(^ose  +  ^^s\n^6.  (41) 

S  =  TT  is  a  critical  point  of  this  function  for  all  values  of  the  parameters,  and  for  gl  <  it  is 

a  local  minimum.  This  may  be  compared  with  the  analysis  of  Example  7,  and  it  is  re-emphasized 
that  the  result  is  independent  of  the  frequency  parameter  A. 


On  the  other  hand,  we  may  explicitly  write  the  fundamental  matrix  solution  $a(^>0)  to  (38)  corre¬ 
sponding  to  the  input  v\.  Assuming  again  that  7  =  ^5,  we  evaluate  this  at  the  fundamental  period 
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t  =  — • 

2/i  f  i-SH.  Ik  \  1 

V  3A€2  a  / 

According  to  Proposition  5,  there  is  a  corresponding  invariant  quadratic  form  specified  by  the 
symmetric  2x2  matrix 

-^V+o(b 

-yf  4h  Jl  +  °<F 

It  is  easy  to  check  that  this  is  positive  definite  for  all  sufficiently  large  A  precisely  when  >  12gi. 
We  summarize  this  observation  in  the  following. 


Proposition  6  Consider  the  cart-pendulum  system  (34)  with  inclination  parameter  value  q  =  | 
and  periodic  forcing  u(')  defined  as  in  (fO)  with  ")  =  fi.  Suppose  the  dynamics  are  linearized  about 
the  equilibrium  6  =  t  of  the  averaged  system  to  yield  (38).  Then  the  origin  is  stable  (in  the  sense 
of  Lyapunov)  for  (38)  under  all  forcing  v\{t)  =  v{Xt)  with  X  sufficiently  large  if  and  only  if  0  =  tt 
is  a  strict  local  minimum  of  the  averaged  potential  (fl)- 


Remark  10  One  might  carry  out  a  similar  stability  analysis  with  the  assumption  7  =  /?  relaxed. 
Any  persistent  bias  of  this  type  in  the  forcing,  however,  seems  to  destroy  the  stability  we  have 
observed.  Thus,  the  shape  of  the  wave-form  forcing  the  system  seems  to  be  crucial.  It  is  worth 
pointing  out  that  numerical  experiments  indicate  that  the  conclusions  of  Proposition  6  remain  valid 
for  any  sufficiently  regular  periodic  input  of  the  form  (24)  (e.g.  vx{t)  —  sinAf). 


Remark  11  Proposition  6  represents  a  remarkable  confluence  of  an  analytical  and  geometrical 
analysis  of  stability.  The  fact  that  local  minima  of  the  averaged  potential  can  be  shown  to  be 
Lyapunov  stable  for  motions  of  the  forced  system  serves  to  further  inspire  the  efforts  we  are  pursuing 
to  simply  characterize  the  nonlinear  response  of  mechanical  systems  to  regular  forms  of  forcing.  It 
is  worth  noting  that  even  for  cases  where  cv  ^  |  and  ^0  is  a  local  minimum  of  the  averaged  potential 
Va  which  does  not  correspond  to  an  eqilibrium  of  the  system  (33),  an  extension  of  the  argument 
used  to  prove  Proposition  6  shows  that  the  linearization  (37)  of  (33)  about  ^0  will  execute  bounded 
motions.  A  more  complete  treatment  of  stability  of  systems  of  the  form  (36)  will  appear  elsewhere. 
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Abstract:  A  large  body  of  recent  literature  has  been  devoted  to  the  topic  of  motions  of  mechanical  systems 
forced  by  oscillatory  inputs.  (See,  e.g.  [7],  [12],  [15],  [17],  and  [29]!)  A  common  feature  in  most  of  this  work 
(with  [7]  being  the  exception)  is  that  the  results  apply  principally  to  kinematic  problems,  and  the  equations  of 
motion  do  not  involve  drift  terms.  The  main  object  of  study  in  the  present  paper  is  the  class  of  bilinear  control 
systems  with  oscillatory  (periodic)  inputs.  We  shall  show  that  included  in  this  class  are  models  of  mechanical 
system  dynamics  which  are  obtained  through  a  process  of  reduction  and  linearization  about  operating  points. 
In  particular,  we  study  the  stability  of  such  systems  under  high  frequency,  periodic  forcing.  The  averaged 
poieniial  is  defined  for  linearizations  of  the  systems  on  the  (reduced)  configuration  space,  and  it  is  shown  that 
stable  motions  of  the  forced  system  are  associated  with  minimum  values  of  this  quantity.  We  use  both  classical 
averaging  theory  as  well  as  a  novel  geometric  argument  to  provide  parallel  but  independent  assessments  of  the 
use  of  the  averaged  poieniial  in  carrying  out  a  stability  analysis.  A  salient  feature  of  the  geometric  approach  is 
that  we  are  able  to  justify  the  use  of  an  energy  like  quantity  to  determine  Lyapunov  stability  in  conservative 
mechanical  systems.  This  makes  contact  with  a  growing  body  of  literature  on  the  use  of  energy  methods  for 
stability  and  control  design  (e.g.  [7],  [8],  [18],  [24],  [26],  [27]).  We  briefly  describe  the  connection  with  previous 
work  on  the  averaged  potential. 


1.  Introduction 


Although  geometric  nonlinear  control  theory  has  achieved,  a  high  degree  of  maturity  over  the  past  two 
decades,  no  substantial  body  of  research  on  nonlinear  control  designs  systematically  utilizing  the  geometry 
of  noncommuting  vector  fields  had  been  reported  until  quite  recently.  (We  cite  Haynes  and  Hermes,  [13], 
as  a  notable  early  effort  along  these  lines,  however.)  Recently  a  number  of  researchers  (e.g.  [12],  [15], 
[17],  [20]  [16],  and  [29])  have  reported  progress  in  the  use  of  open-loop,  oscillatory  controls  to  steer  certain 
classes  of  systems  that  are  of  interest  in  path-planning  problems  in  robotics.  The  present  paper  is  written 
in  the  same  spirit,  but  our  focus  is  on  systems  arising  in  dynamic  rather  than  kinematic  models. 


In  previously  reported  work  ([4]- [7]),  we  studied  stability  of  the  dynamics  of  mechanical  systems  subject 
to  oscillatory  forcing.  The  overall  goal  of  the  research  has  been  to  understand  system  dynamics  which  are 
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the  net  result  of  applying  periodically  varying  forcing  to  a  set  of  noncommuting  vectorfields.  One  result 
has  been  our  characterization  of  stable  operating  points  as  critical  points  of  an  energy-like  quantity  called 
the  averaged  potential  In  particular,  in  [7],  it  was  shown  that  critical  points  of  the  averaged  potential 
which  are  also  strict  local  minima  are  stable  provided  the  frequency  of  the  periodic  forcing  is  sufficiently 
high.  The  results  were  established  using  classical  averaging  theory  together  with  the  adroit  introduction 
of  damping  into  the  models.  In  the  present  paper,  we  show  that  the  stability  of  motions  in  conservative 
dynamical  models  may  be  assessed  in  terms  of  the  averaged  potential  with  no  assumption  of  dissipation  in 
the  system. 


The  paper  is  organized  as  follows.  In  Section  2,  we  study  general  bilinear  systems  forced  by  oscillatory 
inputs.  It  is  shown  that  if  periodic  inputs  are  scaled  by  the  same  factor  in  frequency  and  amplitude,  that 
given  a  certain  assumption  on  average  dissipation  in  the  model,  the  long  term  behavior  will  be  described 
by  an  averaged  model  involving  mean  and  r.m.s.  values  of  the  input  signal.  This  is  a  more-or-less  direct 
consequence  of  the  classical  averaging  theorem.  In  Section  3,  a  special  class  of  mechanical  systems  models 
is  considered.  The  models  are  obtained  through  a  process  of  (Lagrangian)  reduction  and  linearization 
about  operating  points.  Such  models,  with  oscillatory  coefficients,  have  been  treated  extensively  in  the 
Russian  literature.  (See,  for  instance  [11]  and  [30].)  For  the  models  which  are  of  interest  to  us,  it  is  shown 
that  stable  motions  may  be  found  in  terms  of  a  quantity  we  have  called  the  averaged  potential  Section  4 
employs  Floquet  theory  to  prove  that  the  averaged  potential  may  be  used  to  identify  stable  motions  under 
oscillatory  forcing  even  when  dissipation  is  not  present  in  the  model. 


The  use  of  energy  methods  for  analyzing  the  stability  of  mechanical  systems  has  been  pursued  recently  in 
a  number  of  problem  areas  arising  in  the  mechanics  of  rotating  systems  (for  example,  [3],  [8],  [18],  [24], 
[26],  and  [27]).  In  [8]  and  [24],  energy  methods  are  employed  to  design  control  laws  for  the  stabilization  of 
rotating  rigid  bodies.  While  our  idea  of  using  the  averaged  potential  for  studying  the  stability  of  systems 
with  inputs  is  fairly  recent  (See  [4],  [5],  and  [7]  for  details  on  the  averaged  potential)^  there  are  antecedents  of 
the  concept  in  the  analytical  mechanics  literature  from  the  nineteenth  century.  Omohundro  ([22])  describes 
efforts  by  Helmholz,  Hertz,  and  others  to  develop  a  “kinetic  theory  of  matter”  wherein  potential  forces 
would  arise  from  “hidden  or  forgotten  degrees  of  freedom.”  Omohundro  interprets  forces  attributable  to 
such  motions  in  terms  of  quantities  he  calls  pseudo-potentials.  Percival  and  Richards  ([23])  briefly  mention 
such  quantities  under  the  name  effective  potentials. 


2.  Bilinear  Systems 


Consider  the  bilinear  control  system 

m 

X  =  {A  +  '^Ui{t)Bi)x,  (1) 

Z  =  1 

where  A,  5i, . . . ,  Bm  are  constant  n  x  n  matrices,  x{t)  G  R”".  For  much  of  the  analysis  which  follows,  we 
shall  only  assume  that  the  U{(*)’s  are  piecewise  continuous  functions  defined  on  some  interval  [0,//)  where 
0  <  //  <  00.  This  serves  to  distinguish  our  results  from  some  of  the  published  work  on  classical  averaging 
where  it  is  assumed  that  the  forcing  is  continuous  ([25])  or  even  ([14]).  We  shall  also  assume  a  strong 
nilpotency  condition,  namely  that  BiBj  =  0  for  all  z,j  =  l,...,m.  While  this  might  appear  to  limit  the 
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generality  of  our  results,  we  shall  find  that  the  condition  is  always  satisfied  for  a  wide  class  of  systems 
arising  from  models  of  mechanical  system  dynamics. 


Certain  elements  of  the  control  theory  of  such  systems  are  well  developed  (see,  e.g.  [21]),  but  there  are 
few  general  methods  for  design  of  control  laws  to  achieve  particular  motion  objectives.  Exceptions  to  this 
include  the  optimal  control  laws  studied  in  [2]  and  also  the  use  of  oscillatory  inputs  to  approximate  desired 
motions  in  [15]  and  [29].  While  systems  such  as  (1)  are  not  amenable  to  treatment  using  the  methods 
of  [15]  and  [29]  due  to  the  presence  of  the  “drift”  term  A,  these  results  nevertheless  suggest  studying  the 
effect  of  applying  oscillatory  inputs  to  (1). 


Assume  •  •  -  i  lim(')  are  periodic  functions  with  (common)  fundamental  period  T  >  0.  In  order  to  apply 

classical  averaging  theory  to  study  the  motion  of  (1)  we  consider  the  effect  of  increasing  the  frequency  of 
the  forcing.  Specifically,  we  study  the  dynamics  of 

x{t)  =  (A  +  ^  Ui{ijjt)Bi)x{t) 

as  u)  becomes  large.  The  analysis  proceeds  by  scaling  time  and  considering  r  =  ut.  Let  z{t )  =  x{t).  This 
satisfies  the  differential  equation 

^  =  -(A  +  Ui{T)Bi)z.  (2) 

dr  u;  ^ 

Assuming  u  >  0  is  large,  and  writing  €  -  we  see  that  (2)  is  in  a  form  to  which  classical  averaging  theory 
applies. 


Proposition  1  Consider  the  bilinear  system  (1)  with  Ui(-)  —  Ui{^)  where  for  z  —  1, . .  .,m,  zij(')  is  contin¬ 
uous  on  0<t<tj<oo  and  periodic  of  period  T  tj.  Let 


Ui  = 


dt, 


and  let  y{r)  be  a  solution  of  the  constant  coefficient  linear  system 

y{T)  =  €(A  +  ^  UiBi)y{T).  (3) 

If  zq  and  yo  are  respective  initial  conditions  associated  with  (2)  and  (8)  such  that  \zo  —  yo\  —  ^(f)>  then 
\z{t)  —  y(r)|  =  0(e)  on  a  time  scale  r  ~ 


This  proposition  is  a  straightforward  application  of  classical  averaging  theory;  see  [25].  Unfortunately,  it 
doesn’t  capture  the  full  range  of  behaviors  that  are  possible  with  different  scalings  of  the  magnitudes  of 
the  inputs  ii,-.  For  instance,  consider  the  following  second  order  system 

£•(<)  + (a  + /3u(t))a:(t)  =  0.  (4) 

When  written  in  first  order  form,  this  is  a  special  case  of  systems  of  the  form  (1).  If  we  simultaneously 
scale  the  frequency  and  magnitude,  u(t)  uu(ut),  the  first  order  form  of  (4)  becomes 

/  0  lWxi\ 

\X2  )  \—a  —  f3uu(ujt)  0y\X2/’ 


(5) 
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Figure  1:  Stability  of  Hill’s  equation  (4)  with  square-wave  forcing  (7). 


and  direct  appeal  to  Proposition  1  is  not  possible.  Nevertheless,  re-writing  the  second  order  system  in 
terms  of  the  “slower”  time  scale  r  —  ut  and  z{r)  =  x{t),  we  obtain 


d?z 


a  (3 
—  H — u 
u 


=  0. 


(6) 


When  ri(r)  cost,  (6)  is  a  standard  rendering  of  Mathieu’s  equation  (See,  e.g.  [28].),  and  in  this  case,  the 
classical  theory  states  that  the  origin  is  Lyapunov  stable  for  sufficiently  large  uj  provided  a  >  —  This 
result  will  appear  as  a  special  case  of  the  theory  to  be  developed  in  the  remainder  of  the  paper. 


Remark  1  The  stability  of  (4)  under  periodic  forcing  has  been  most  widely  studied  for  the  specific  case 
u{t)  =  cost.  There  is  growing  interest  in  other  types  of  inputs.  (See,  e.g.,  [19].)  For  square  wave  input 


w(0  = 


-1,  if  0  <  <  <  ^; 
1,  if  ^  <  <  <  1, 


(7) 


the  regions  in  the  (a,/3)  parameter  space  for  which  the  origin  (x,a;)  =  (0,0)  is  a  stable  rest  point  for 
(4)  are  depicted  in  Figure  1.  Qualitatively,  this  picture  is  quite  similar  to  the  classical  rendering  of  the 
stability  regions  corresponding  to  the  forcing  u{t)  =  cost.  The  boundary  of  the  left-most  stability  region 
is  approximately  parabolic  near  the  origin  and  satisfies  a  =  +  ^(/3^)*  (CL  [28],  pp.  208-213.) 


A  particularly  simple  approach  to  averaging  systems  with  high-frequency  and  high-amplitude  inputs  in¬ 
volves  making  a  time-varying  change  of  coordinates  which  replaces  the  inputs  Ui{-)  with  their  integrals. 
The  coordinate  transformation  will  turn  out  to  have  an  especially  simple  form  in  the  case  that  BiBj  =  0 
for  aU  pairs  i,^j  =  1, . .  .,m.  Moreover,  this  condition  will  always  be  satisfied  by  the  mechanical  systems  to 
be  considered  in  the  next  section.  We  shall  now  further  suppose  that  each  Ut(*)  is  a  zero-mean  periodic 
function  of  period  T.  This  assumption  involves  no  loss  of  generality.  If  we  wish  to  study  a  system  (1)  in 
which  the  mean  of  the  inputs  i/t(*)  is  not  zero,  we  may  write  Ui{t)  =  Ui  +  where  Ui  —  ^  Jq  ds. 
The  dynamics  of  (1)  will  be  the  same  as  for  the  system  x  =  {A  +  where  A  =  A  +  YluiBi  and 

which  satisfies  our  assumption  on  the  mean  of  the  inputs.  Letting 
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it  is  easy  to  see  that  each  Vi(-)  is  also  a  periodic  function  of  period  T.  Given  that  a;(-)  satisfies  (1),  we 
transform  to  new  coordinates  q(-)  using  the  formula 

9(0  =  (^  -  E  Vi{t)Bi)x{t).  (8) 

A  straightforward  calculation  shows  that 

q{t)  =  (a  +  '^Vi{t)[A,Bi]  -Y,'>’iit)vj{t)BiABj)q{t),  (9) 

where  [A^Bi]  denotes  the  usual  matrix  Lie  bracket  ABi  —  BiA. 


The  utility  of  changing  to  ^-coordinates  follows  from  the  fact  that  in  the  slow  time  scale  r  =  ut,  (9) 
assumes  a  form  to  which  classical  averaging  theory  applies,  provided  the  scaling  factor  cj  is  taken  to  be 
sufficiently  large. 


Theorem  1  Consider  the  dynamical  system 


m 

i  +  y^u!Ui{ut)Bi^x, 


2  =  1 


(10) 


where  for  i  —  1, . . . ,  m  Ui{*)  is  a  piecewise  continuous  periodic  function  with  period  T  and  mean,  0,  arid  for 
j  =  1, . . m,  BiBj  —  0.  Let  c  =  Define  for  each  i  =  1, . . . ,  m,  periodic  function  Vi{t)  =  Jq  Ui(s)  ds, 


and  let 


Let  y{t)  be  a  solution  of  the  constant  coefficient  linear  system 

y  =  (^  +  -  <^ij)BiABj'jy.  (11) 

hi 


Suppose  the  eigenvalues  of  (11)  have  negative  real  parts.  Then  there  is  a  ti  >  0  such  that  for  u  >  0 
sufficiently  large  (i.e.  for  e  >  0  sufficiently  small),  if  xo  and  yo  are  respective  initial  conditions  for  (10) 
and  (11)  such  that  |xo  —  j/o|  =  0{e),  then  |x(t)  —  y{t)\  =  0{e)  for  all  t  >  t^. 


Proof:  Given  the  system  (10)  and  the  functions  Uj(-)  defined  in  the  statement  of  the  theorem,  define  the 
new  state  coordinates 

m 

9(<)  =  a  - 

2  =  1 

in  terms  of  which  the  dynamics  are  expressed  by  (9).  Changing  to  time  variable  t  =  uit  and  state  variable 
z{t)  =  q{t),  we  have  dynamics 


^  =  e(A  +  5]t;,(r)[A,fii]  -  J]ui(r)uj(r)HiABj)^. 


(12) 
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It  is  this  system  to  which  the  averaging  theorem  ([25],  Theorem  4.2.1,  p.  71)  applies,  with  the  corresponding 
averaged  dynamics  given  by 


dz 

dr 


e(A  +  ^  Vi[A,  Bi]  -  <yijBiABj^  z. 


(13) 


Specifically,  if  zq,  zq  are  respective  initial  conditions  of  (12)  and  (13)  such  that  \zo  —  fo|  =  0(e),  then 
\z(t)  —  z(t)\  —  0(e),  on  a  time  scale  r  ~  If  the  spectrum  of  (13)  is  in  the  left  half-plane,  then  the 
conclusion  is  valid  on  the  time  interval  [0,oo). 


Reverting  to  the  original  time-scale,  and  writing  q(t)  =  z(ut),  we  find  that  q  satisfies 

9  =  (^  +  S  ~  S  ^ijBiABj^  q(t), 

and  given  \q{0)  —  ^(0)|  =  it  follows  that  \q{t)  -  ^(/)|  =  0(e),  and  provided  the  spectrum  is  in  the  left 

half-plane,  this  is  valid  on  the  time  scale  [0,oo). 


Let  x(t)  =  (/ +  J2viBi)q(t),  This  satisfies  the  constant  coefficient  linear  differential  equation 

X  ~  (^ A  +  ^^(viVj  ~  crij)BiABj^x, 
ij 

Noting  also  that  x(t)  =  (I  +  2]  we  write 

kV)-^(0l  =  \(I  +  ^ViBi)q(t)-(I  +  J2vrB,)q(t)\ 

<  |9(t)  -  q(t)\  +  1  ViB,q(t)  -  Y,  ViBiq(t)  +  ~  ^i)Biq(t)\ 

<  |(/(t)  -  ^(<)|  4-  \YviBi(q(t)  -  9(0)1  +  “  ^i\\Biq(t)\ 

<  -^i|9(0  -  9(01  +  ■^2|9(0I> 

where  Mi  =  supo<^<j(l  +  l'yf(/)5i|)  and  M2  =  ^^Po<t<T(l^  \vi{t)-Vi\\Bi\).  The  norms  of  vector  quantities 
in  this  case  are  the  ordinary  finite  dimensional  Euclidean  norms,  and  the  matrix  norms  are  the  induced 
norms.  Under  the  hypothesis  that  the  spectrum  of  (11)  lies  in  the  left  half-plane,  it  is  also  the  case  that  the 
spectrum  of  (13)  lies  in  the  left  half-plane.  Hence,  from  the  above,  |?(f)  g(i)|  =  0(e),  for  c  >  0  sufficiently 

small.  Moreover,  if  the  eigenvalues  of  (13)  are  in  the  left  half-plane,  it  follows  that  the  equilibrium  point 
q  =  0  IS  asymptotically  stable,  and  |g(/)|  ^  0  as  ^  increases.  This  proves  the  theorem.  □ 


Remark  2  Although  Theorem  1  is  of  some  interest  in  characterizing  the  response  of  bilinear  systems 
to  periodic  forcing,  a  limitation  may  be  uncovered  by  re-examining  (5).  Again  suppose  u(t)  =  cost. 
Then  T  =  27r  and  v(t)  =  sinL  Hence  v  =  0,  a  =  1/2,  and  the  coefficient  matrix  appearing  in  (11)  is 

ip2  *  Since  the  eigenvalues  are  purely  imaginary.  Theorem  1  does  not  apply.  On  the  other 
hand,  as  remarked  earlier,  the  origin  is  stable  in  the  sense  of  Lyapunov  for  all  sufficiently  large  values  of  u. 
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3.  Linear  Lagrangian  Systems 


Consider  a  Lagrangian  control  system  prescribed  by  a  Lagrangian 

1  ^  m 

L{x,  v)  =  -x'^Mx  +  x^Sx  -  -x'^Kx  +  ^  ViX^CiX 

2  ^  i^\ 

where  M,  K,  and  Ci  are  n  X  tz  symmetric  matrices  with  M  and  K  positive  definite,  and  5  is  an  n  x  n 
skew-symmetric  matrix.  (See  [21]  and  [7]  for  details  on  Lagrangian  control  systems.)  The  corresponding 
dynamics  are 

Mx  +  2Sx  +  Kx  +  =  0  (^4) 

where  Ui{i)  =  Vi{i)  for  i  =  1, . . . ,  m. 


Remark  3  Linear  Lagrangian  systems  arising  from  a  class  of  super- articulated  mechanical  system  models. 
Linear  Lagrangian  control  systems  of  the  form  (14)  arise  as  partially  linearized  models  of  so-called  super- 
articulated  mechduuc^l  systems,  (See  [4]-[7]  for  an  introduction  to  super-articulated  mechanical  systems.) 
The  starting  point  for  a  study  is  a  Lagrangian  L[q^  q)  defined  on  the  tangent  bundle  TQ  of  the  configuration 
space  of  a  mechanical  system.  Suppose  that  the  generalized  coordinates  can  be  partitioned  3.s  q  =  (<jfi?^2) 
and  that  exogenous  generalized  forces  (control  inputs)  can  be  applied  to  only  the  coordinates  ^2,  while  the 
coordinate  variables  comprising  qi  evolve  freely.  Specifically,  we  assume  the  equations  of  motion  for  the 
system  take  the  form 

(15) 

dt  dqi  dqi 


d  dL  dL 
dt  dq2  dq2 


(16) 


where  z/  is  a  p-vector  of  controls,  and  we  suppose  the  mapping  92  ^  is  invertible.  For  the  purposes 
of  the  model,  we  assume  that  the  components  Ui{-)  are  each  piecewise  analytic  functions  on  [0,oo).  For  a 
wide  variety  of  interesting  physical  systems,  the  terms  on  the  left  hand  side  of  (15)  depend  explicitly  on 
the  generalized  velocity  52  but  not  on  the  generalized  coordinate  q2  itself.  It  is  to  such  models  that  we 
direct  our  attention. 


The  problem  we  wish  to  consider  is  that  of  specifying  u(-)  to  steer  the  state  (qi^qi)  in  some  desired 
fashion.  While  it  is  natural  to  view  the  velocities  ^2  ^.s  intermediate  variables  which  transmit  the  effects  of 
the  controls  to  the  states  via  the  coupled  equations  (15)- (16),  it  will  more  be  useful  in  what  follows  to  view 
the  variables  ^2  themselves  as  inputs,  with  the  states  (^1,^1)  being  determined  by  (15)  and  the  required 
generalized  forces  t/(-)  determined  by  (16). 


To  treat  a  fairly  general  class  of  systems  which  will  lead  to  models  of  the  prescribed  type,  we  further 
partition  52  =  (921? 922)5  assume  that  the  velocities  ^21  ^re  held  to  certain  constant  values  while 
922  =  '^^(0  is  assumed  to  be  a  vector  of  inputs,  trajectories  of  which  may  be  chosen  to  control  the  system. 
By  holding  certain  velocities  (921)  constant,  it  will  be  seen  that  gyroscopic  coupling  is  introduced  among 
the  components  of  the  states  (91,91).  Physical  examples  of  such  coupling  may  be  found  in  multibody 
system  in  which  certain  component  bodies  are  driven  to  turn  at  a  constant  angular  velocity  relative  to 
some  body  coordinate  system.  The  mechanics  of  systems  with  gyroscopic  coupling  have  been  studied  by 
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a  number  of  researchers.  See,  for  example  Bloch,  et  al.  [9].  For  simplicity  of  exposition,  we  shall  assume 

921  =  1  =  ,1)^. 


Suppose  the  Lagrangian  has  the  form  L{q,q)  =  \q^M{q)q  -  V{q)  where  M{q)  may  be  partitioned  con¬ 
formably  with  (91^,921,922)  3-^ 

/  Mii{qi)  Mi2(9i)  Mi3(qi)\ 

M{q)  =  Mi^(9i)  M22iqi)  M23 

^^23  -^^33  / 

We  assume  91  is  n-dimensional,  as  is  x  in  (14),  921  is  ^'-dimensional  for  some  ly  >  0,  and  922  is  m-dimensional, 
the  dimension  of  the  control  input  space  in  (14).  Then  (15)  may  be  rendered 

+  (•■^0(91)  -  Jo(qi)^)4i  +  +  ^k(i)(^k(qi)  -  Jk(qi)^)qi]  =  0  (17) 

where  i/?  is  the  reduced  Lagrangian 

1  r 

^h(9i,9i)  =  2?!  -^^11(91)91  -  yR{qi) 

with  Fh(9i)  =  l'^(9i)  —  ^1^M22(9i)1.  Jo{qi)  is  the  n  x  n  matrix  whose  (i,j)-th  entry  is  the  partial 
derivative  with  respect  to  qj  of  the  sum  of  entries  in  the  i-th  row  of  the  sub-block  Mi2(9i),  Mj^3{qi)  is  the 

^'-th  column  of  the  sub-block  Mi3{qi),  and  Jkiqi)  is  the  Jacobian  matrix  .  Suppose  910  €  R”  and 

... 

=  (910? 0)  is  a  rest  point  of  (15)  for  all  values  of  q2.  If  we  expand  the  dynamics  (17)  about  the  rest 
point  and  retain  terms  which  are  jointly  first  order  in  we  obtain 

M9  -f  59  -I-  ii:'9  -I-  Yl['^k{t)Ckq  +  Vk{Ck  -  Cl)q]  =  0,  (18) 


where  9  =  ^91  is  the  variational  term,  M  =  Mii(9io),  K  =  ^^(910),  5  =  Jo(9io)  -  <^o^(9io),  and 

Cfc  =  ^(910). 


Note  that  in  (18),  the  matrices  are  not  necessarily  symmetric.  The  presence  of  terms  Ck  —  C'^  ^  0 
indicates  a  form  of  “inertial  coupling”  between  the  inputs  v  and  the  states  which  is  not  present  in  (14), 
and  systems  of  the  form  (18)  thus  appear  to  be  more  general  and  to  include  (14)  as  a  special  case.  Indeed, 
the  case  in  which  the  above  linearization  leads  to  symmetric  coefficient  matrices  Ci  is  somewhat  special, 
although  it  includes  a  number  of  interesting  sub-cases  including  systems  in  which  dim  qi  —  1,  There 
appear  to  be  some  fundamental  differences  between  the  cases  of  symmetric  and  non-symmetric  C^’s.  While 
a  parallel  theory  of  systems  with  non-symmetric  Cj’s  as  described  in  this  section  (where  the  models  include 
the  effects  of  dissipation)  could  be  described,  we  shall  not  do  this.  The  results  of  the  following  section  do 
not  generalize  as  easily,  and  we  shall  conclude  the  paper  with  a  remark  detailing  the  differences  between 
the  two  cases.  □ 


Example  1  A  well  known  physical  phenomenon  which  may  be  used  to  illustrate  our  theory  is  stabilization 
of  the  inverted  pendulum  by  means  of  oscillatory  forcing  of  the  point  of  suspension.  (Cf.  Arnold  et  al.,  [1], 
pp.  152-153.)  The  equation  for  this  system  is 
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where  t  is  the  length,  g  is  the  gravitational  acceleration,  and  y{t)  is  the  vertical  height  of  the  point  of 
suspension  at  time  t.  While  y  could  be  viewed  as  a  physical  coordinate,  we  shall  instead  view  i/  as  a 
control  input,  which  we  rename  u.  (9,6)  =  (7r,0)  is  a  rest  point  for  the  system  for  all  choices  of  u,  and 
carrying  out  the  above  linearization,  (14)  specializes  to 

=  (19) 


□ 


To  make  contact  with  the  preceding  section,  we  rewrite  equation  (14)  in  first-order  form: 


0  I 

-M-^K 


;) +i;  «.(<)(_, 1) 


(20) 


where  Xi  =  x  and  X2  =  i.  One  would  expect  that  systems  of  this  form  do  not  meet  the  hypotheses  of 
Theorem  1,  since  such  models  apply  to  conservative  mechanical  systems.  It  is  fairly  standard,  however, 
to  introduce  dissipation  to  such  models  by  means  of  a  quadratic  (Rayleigh)  dissipation  function  ^x^ Dx, 
where  L)  is  a  positive  definite  symmetric  matrix.  With  such  dissipation  included  (20)  is  rewritten  as 


/  0  I  ^ 

-M-\2S  +  D))  \X2 


i-l 


0  0\  fxi\ 

-M-^Ci  q)\X2)' 


(21) 


For  i  =  1,. .  .,m,  let  Wi{-)  be  a  piecewise  continuous  periodic  function  of  period  T,  having  zero  mean,  and 
consider  the  dynamical  system  (20)  (resp.  (21))  forced  by  u,(/)  =  u’Wi{ujt)  for  some  u  >  0.  Both  (20)  and 
(21)  have  the  form  (10).  The  averaged  system  corresponding  to  (21)  takes  the  form 

(I)  "  ( -M-i A  -h  Ei,j(viVj  -  cTij)M-'^CiM-^Cj  -M-\2S  +  D))  {y 2)  ' 

(Cf.  (11).)  Note  that  this  system  is  associated  with  the  constant  coefficient  Lagrangian 

Hy,y)  =  ^y^My  -  ~  ViVj)CiM~'^Cj^y  +  fSy, 

with  Rayleigh  dissipation  function  ^y^Dy.  The  quantity 

i/(A  +  -  ViVj)CiM-^Cj)y 

plays  the  role  of  a  potential  and  we  are  lead  to  the  following. 


Definition  1  For  the  Lagrangian  system  (21)  with  periodic  inputs  Ui{t)  =  uwi(ut)  as  prescribed  above 
and  given  the  related  averaged  quantities  Vi  {i  =  1, . . . ,  m)  and  aij  {i,j  =  1, . . . ,  m)  as  defined  in  Theorem 
1,  we  define  the  averaged  potential 

Va{x)  =  ix^(A  -1-  -  ViVj)CiM-^C^)x. 
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Remark  4  It  is  easy  to  see  that  this  quantity  coincides  with  the  averaged  potential  defined  previously  for 
mechanical  systems  in  [7] 


As  in  [7],  the  averaged  potential  may  be  used  to  determine  the  stability  characteristics  of  the  dynamical 
system  (21)  under  the  influence  of  oscillatory  forcing.  The  following  theorem  gives  the  result  more  precisely. 


Theorem  2  Let  be  a  piecewise  continuous  periodic  function  of  period  T  and  mean  0  for  i  ~  1, . . . ,  m. 
Define  for  each  i  =  1, . . m,  the  periodic  function  Vi{t)  ~  ds,  and  let 


1 

Vi  =  -J^  Vi{s)(ls, 


i  =  1 , . . . ,  ??i 


1  /-^ 

Vi{s)vj{s)  ds,  tj  =1,...,  m. 


Suppose  D  is  positive  definitej  and  let  the  system  (21)  be  forced  by  the  oscillatory  inputs  Ui{t)  =  LJWi{ujt) 
for  f  =  1, . . . ,  m.  If  the  averaged  potential 


-  vrvj)CiM-^C,): 


is  a  positive  definite  quadratic  form,  the  state  (.x'i,a;2)  =  (0?0)  (asymptotically)  stable  for  all  sufficiently 
large  oj. 


Proof:  Suppose  the  averaged  potential  is  positive  definite.  Then  the  “energy” 

E  =  \ylMy2  +  \yj {k  +  -  ViVj)CiM-^Cj)yi 

is  a  positive  definite  quadratic  form  on  (t/i,  2/2)-space.  Differentiating  this  along  trajectories  of  (22),  we 
obtain 

Thus  E  serves  as  a  Lyapunov  function,  and  invoking  LaSalle’s  invariance  principle,  all  trajectories  of  (22) 
approach  the  set  of  points  where  ?/2  =  0-  Since  the  only  solution  to  (22)  with  t/2  =  0  is  (j/i,  2/2)  =  (0»0),  we 
conclude  that  the  origin  is  an  asymptotically  stable  equilibrium.  Hence  the  eigenvalues  of  (22)  must  lie  in 
the  left  half-plane,  and  the  conclusion  of  the  theorem  follows  from  Theorem  1.  □ 

Theorem  2  relies  on  classical  averaging  theory  to  establish  the  connection  between  the  averaged  potential 
and  stability  of  forced  systems.  It  will  be  shown  in  the  next  section  that  the  result  is  geometric  in  nature, 
and  the  stability  analysis  may  be  applied  to  conservative  systems  as  well. 
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4.  A  Geometric  Theory  of  Averaging 


As  pointed  out  in  [7],  the  averaged  potential  provides  useful  information  about  the  stability  of  motion 
even  in  the  absence  of  dissipation.  A  natural  tool  in  this  study  is  Floquet  theory,  and  while  it  may  be 
inconvenient  to  analyze  an  explicit  representation  of  the  time-T  mapping  associated  with  a  system  of  the 
form  (10),  the  following  lemma  provides  a  useful  approximation.  The  hypotheses  adopted  in  this  section 
are  essentially  the  same  as  in  the  preceding  two,  but  no  assumptions  are  made  regarding  the  asymptotic 
stabihty  of  the  origin  for  the  averaged  system. 


Lemma  1  Consider  the  dynamical  system 


i  =  (a  +  '^u>Wi{ujt)Bi^x, 

t=i 

where  for  i  =  1,. .  .,m  Wi^r)  is  a  piecewise  continuous  periodic  function  with  period  T  and  mean  0  and 
continuous  from  the  right  for  all  r  €  [0,T].  Suppose  further  that  for  i,j  =  1, . . . ,  m,  BiBj  =  0.  Let  e  =  K 
Define  for  each  i  —  1, . .  .,m,  the  periodic  function  n,(t)  =  u;,(s)  ds.  Then 


#4T,  0)  =  /  +  ^  exp  [-  ^  Vi{g)Bi^  •  A  •  exp  Vi{r])Bi  dp  •  e  +  o{e) 
is  the  transition  matrix  of 


Y~  +  X!  MT)Bi)  y{T). 


dr 


(23) 


Proof:  Let  t  =  ut  and  define  ?/(r)  ==  x{t).  If  $£(r,  s)  is  the  transition  matrix  of  the  system  (23),  then 

(Ae  +  ^u;.(r)R.)$£(r,s), 

with  the  derivative  being  taken  from  the  right  at  points  of  discontinuity  of  any  of  the  Differentiating 
both  sides  with  respect  to  c,  we  find 

IJi  =  A$e(r,  s)  +  (Ac  +  ^  Wi{T)Bi)  ^(r,  s). 

It  follows  from  the  variation  of  constants  formula  that 

^(r,  s)  =  $£(r,  s)^(5,  s)  +  ^  V)  A  ^e(v, «)  dp. 

Since  s)  =  0,  the  first  term  on  the  right-hand  side  of  this  equation  vanishes. 

Letting  e  =  0,  we  see  that  $o('r,s)  is  just  the  transition  matrix  associated  with  the  coefficient  matrix 
and  under  the  assumption  BiBj  =  0  for  all  f,  j,  we  may  write 

$o(r,  s)  =  exp  l/T  Wi(p)Bi  dp\ . 
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Hence 


de  I, 


=  ■  A-exi^y^  '^Wi{^)Bid^^d7]. 


Note  also,  since  we  have  assumed  each  u;,(-)  to  have  zero  mean,  we  have  #0(^5  0)  =  I.  The  formula  for  the 
time-T  map  may  thus  be  given  in  terms  of  a  series  expansion  with  respect  to  c: 


$.(r,o)  = 


This  proves  the  lemma. 


rT  rT  ni) 

^  ^  Jo  f  Jr,  ^  d^j-A-  exp  y  ^  'Wi{^)Bi  dr]  ■  e  +  o(c) 

^  ^  Jo  ^  Mv)Bi\  •  A  •  exp  Vi{T])Bi^  dr]  ■  e  +  o(e). 


□ 


For  linear  Lagrangian  systems  without  damping,  the  expansion  derived  in  Lemma  1  corresponding  to 


0  / 

-M-^K  -2M-15 


+  '^UWi{u!t) 


i=l 


0 


takes  the  form 


(24) 


^(T0)=f  ViM-^CiTc  IT(  \  . 

-  ■2J:v,SM-^Q  +  E  (T,,CiM-^C,)Te  I  -  2M-^S  +  E  m-^CiTe  j  + 

(25) 

The  stability  analysis  of  (24)  will  be  carried  out  in  terms  of  iterates  of  the  time-T  transition  matrix, 
$e(r,0).  More  specifically,  it  follows  from  Floquet  theory  that  the  origin  will  be  stable  under  the  motion 
(24)  when  iterates  of  $e(r, 0)  are  stable.  (See  Guckenheimer  and  Holmes,  [14],  pp.  24-25,  for  a  discussion 
of  the  relevant  Floquet  theory.) 


Lemma  2  Let  $e(T,  0)  be  the  transition  matrix  associated  with  (24)  via  equation  (25)^  where  T  is  the 
common  fundamental  period  of  the  Wi{*)’s.  0)  has  distinct  eigenvalues  lying  on  the  unit  circle^ 

then  there  is  a  symmetric  positive  definite  matrix  M  such  that 

$,(r,0)^Af$,(r,0)  =  A4.  (26) 

Conversely^  if  $e(T,  0)  leaves  a  symmetric  positive  definite  M  invariant  as  in  (26),  then  has 

eigenvalues  on  the  unit  circle,  and  iterates  $e(r,  0)^  remain  hounded  as  n  ^  oo. 


Proof:  If  all  eigenvalues  of  $(r,  0)  are  distinct  and  lie  on  the  unit  circle,  there  is  a  nonsingular  matrix  P 


such  that  A  =  P  ^$(7",  0)P  is  block  diagonal  with  O’s  and  2x2  blocks  of  the  form 


cos  (j)  —  sin  0 


sin  4>  cos  (j) 

the  diagonal.  Now  AA^  =  /,  and  from  this  it  follows  that  $(r,  0)^A^$(T,  0)  =  M  where  =  PP^. 


on 


On  the  other  hand,  if  $e(T, 0)  leaves  invariant  a  symmetric  positive  definite  M,  let  be  the  unique 

positive  definite  symmetric  square  root  of  M.  It  is  easy  to  see  that  A4^/^$c(T,  0)A4“^/^  is  an  orthogonal 
matrix.  Hence  A4^/^$e(r,  has  eigenvalues  on  the  unit  circle,  and  the  same  conclusion  for  $e(r,  0) 

is  an  immediate  consequence.  □ 


The  remaining  results  complete  our  stability  theory  by  establishing  conditions  under  which  $e(r,  0)  given 
in  (25)  leaves  a  symmetric  positive  definite  quadratic  form  invariant. 
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Lemma  3  For  the 
blocks 


system  (24)  evolving  in  R^”,  let  ^^(T,  0)  as  given  in  (25)  be  partitioned  into  four  nxn 


\$21 


$12 

$22/ 


Suppose  $£  is  symplectic,  satisfying  $f  J$£  =  J ,  where  J  = 


0  I 
-I  0 


is  the  standard  symplectic  bilinear 


form  on  R^".  Then  $£  also  leaves  a  symmetric  bilinear  form  M  invariant.  This  form  may  be  explicitly 

written  ,  ^ 

-$2i  -  $21  $11  ~  $22 

$11  —  $^2  $12  +  $12  . 


M  = 


(27) 


Proof;  The  proof  is  a  straightforward  verification  that  the  linear  (in  M)  equation  $  =  Af  is  satisfied. 
□ 


Theorem  3  Consider  (21)  and  the  corresponding  averaged  system  (22)  (with  D  =  Oj.  Assume  that  the 
inputs  Ui{t)  =  u)Wi{ojt)  and  the  averaged  quantities  Vi  and  aij  are  as  in  Theorem  2,  and  suppose  further 
that  5  =  0.  If  the  averaged  potential  is  positive  definite,  the  origin  is  stable  in  the  sense  of  Lyapunov  under 
the  motion  of  (20)  provided  u)  is  sufficiently  large. 


Proof:  Change  variables  to  Zi  =  M^Xi  where  Mz  is  the  unique  symmetric,  positive-definite  square-root 
of  the  matrix  M.  The  equation  (20)  becomes 


^2 


(28) 


where  K  =  AM“2,  Ci  =  Letting  u,(f)  =  uwi{Lvt),  the  time-T  map  of  Lemma  1  for 

this  system  takes  the  form 


$e(T,0)  = 


I  0 
0  I 


+ 


—(A  +  E  E 


Te-f  o(e). 


We  must  find  necessary  and  sufficient  conditions  under  which  it  is  possible  to  to  solve  —  Ai  for  a 

symmetric,  positive  definite  M.  Since  all  coefficient  matrices  appearing  in  (28)  belong  to  the  symplectic 
Lie  algebra  sp{n,K),  $£  =  $c(r,0)  belongs  to  the  group  of  2n  X  2n  symplectic  matrices,  which  satisfy 
=  J.  Hence  from  Lemma  3,  we  conclude  that  $£(T,  0)  also  leaves  the  symmetric  quadratic  form 
(27)  invariant.  Using  this  explicit  representation,  we  write 


M,= 


(  K  -b  Yfc^ijCiCj 
\  -EviCi 


-  E 

I 


Now  Mi  will  be  positive  definite  for  all  sufficiently  small  e  (large  w)  precisely  when  the  the  matrix 

V  -EviCi  I  J 


is  positive  definite.  This  is  more  easily  analyzed  in  the  coordinate  system 
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with  the  transformation  from  ^  to  i  changing  the  expression  of  the  matrix  in  question  to 

^  K  +  -  ViVj)CiCj  0  ^ 

The  positive  definiteness  of  this  matrix  is  clearly  equivalent  to  the  positive  definiteness  of  the  sub-block 
+  E(^»i  -  ViVj)CiCj.  This  matrix  is  obviously  positive  definite  if  and  only  if  the  a;-coordinate  averaged 
potential  matrix  K  -|-  “  ViVj)CiM~'^Cj  is  positive  definite.  This  proves  the  theorem.  □ 


Remark  5  The  assumption  that  6’  =  0  appears  to  be  needed  for  the  construction  of  the  invariant  quadratic 
form  in  terms  of  which  we  determine  the  stability  of  this  system.  Beyond  this,  however,  it  has  been  shown 
that  the  stability  of  linear  systems  with  gyroscopic  terms  is  ambiguous.  In  particular,  it  was  shown  in  [6] 
and  also  [9]  that  in  cases  where  K  (or  K is  not  positive  definite,  the  system  (14) 
will  have  right  half  plane  eigenvalues  if  small  amounts  of  dissipation  are  formally  introduced,  as  was  done 
in  Section  3.  Thus,  the  conditions  of  the  theorem  would  stand  no  chance  of  being  necessary  and  sufficient 
in  the  presence  of  non-zero  S. 


Example  2  (Example  1,  reprise.)  In  Example  1,  the  linear  Lagrangian  system  corresponding  to  the 
inverted  pendulum  with  oscillatory  support,  was  written  (19): 

Identifying  the  standard  terms  in  (14)  for  this  example  we  have  K  =  -g/i,  M  =  and  C  =  —1/i.  For 
the  oscillatory  input  u{t)  ~  cj^sin(cji),  we  have  u  =  0  and  The  averaged  potential  for  this 

system  is  The  criterion  for  stability  is  that  the  point  0  =  0  is  a  strict  local  minimum  of  this 

function,  and  this  will  obviously  be  the  case  precisely  when  >  2g£.  This,  of  course,  coincides  with  the 
classical  stability  result,  (Cf.  [1],  p.  153.)  for  this  problem.  □ 


Example  3  Consider  the  electric  circuit  in  Figure  2.  The  current  in  the  circuit  satisfies  the  second  order 
differential  equation 

L^I{t)  +  u(t)I{t)  =  0 

where  u(t)  depends  on  the  position  of  the  switch.  Specifically,  suppose,  the  switch  obeys  a  periodic 
switching  law 

ri/Ci,  ifO</<A/i; 
u{t)  =  ll/C2,  if  Xh<t<h, 

[  Extend  to  be  h-periodic, 

where  0  <  A  <  1  determine  the  fraction  of  each  duty  cycle  the  switch  is  closed  to  the  left.  Following  the 
methods  of  the  paper,  we  can  write  down  the  averaged  potential  for  this  system.  This  is  reflected  in  the 
evolution  equation  for  the  averaged  system: 


where 


^I+(u/L  +  f3^/L^)!=0, 


T  +  and  =  h^l/C^  -  1/C,)\b 


A2(A-  1)2 


12 
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Figure  2:  A  switched  LC-circuit. 

In  this  case,  the  averaged  potential  together  with  the  mean  of  the  input  reflects  the  effective  capacitance 
of  the  switched  circuit. 


5.  Concluding  Remarks  on  Connections  with  Problems  in  Mechanics 


While  the  theory  we  have  developed  applies  directly  to  provide  an  approach  to  stability  analysis  for  classical 
bilinear  control  systems,  and  thus  applies  to  models  describing  electrical  networks,  it  has  also  been  shown  in 
Section  3,  that  certain  models  arising  Lagrangian  mechanics  can  also  be  analysed  by  our  averaging  method, 
and  in  particular  stability  could  be  analyzed  in  terms  of  an  energy-like  quantity  call  the  averaged  potential. 
In  Section  3,  we  imposed  a  certain  symmetry  restriction  on  the  mechanical  systems  being  considered.  We 
conclude  by  discussing  how  this  restriction  can  be  relaxed. 


Remark  6  (The  case  of  non-symnietric  Ci.):  In  the  case  of  a  system  of  the  form  (18),  in  which  the  Cfs 
are  not  symmetric,  we  may  nevertheless  apply  the  coordinate  transformation  (8),  which  in  this  case  is 
rendered  r.\/\ 

ill)  "  /)(a;2)' 

When  5  =  0,  the  dynamics,  written  in  terms  of  ^-coordinates  evolve  according  to 


-M-.A'  o) 


0 


2)1  C: 


Averaging  the  coefficients  over  one  period  and  applying  the  “inverse”  coordinate  transformation, 

ill)  ^  /)(gS)’ 

the  “averaged  system”  dynamics  take  the  form 

(^2)  "  {-M-\K  -1-  -M-\Evi{Ci  -  Cj)))  (2/2)  ■ 
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The  quadratic  form  defining  the  averaged  potential  for  this  system  is  thus  seen  to  be  K  +  ~~ 

ViVj)Cj' If  this  is  positive  definite,  stability  results  along  the  lines  of  Theorem  2  may  be  ob¬ 
tained  for  the  case  in  which  a  positive  definite  damping  term  (associated,  say,  with  a  Rayleigh  dissipation 
function  as  described  in  Section  3)  enters  the  dynamics.  For  the  case  in  which  there  is  no  dissipation 
present,  the  constructions  presented  in  Section  4  are  less  straightforward,  and  the  stability  theory  for  the 
case  of  non-symmetric  C^’s  remains  incomplete.  □ 


Remark  7  (The  relationship  with  earlier  work  on  the  averaged  potential.)  In  [7],  we  introduced  a  more 
general  version  of  the  averaged  potential  aimed  at  analyzing  the  motions  of  a  large  class  of  classical  me¬ 
chanical  systems.  For  this  larger  class  of  systems,  it  can  be  shown  that  the  results  of  the  present  paper 
may  be  applied  directly  to  linearizations  of  the  dynamics  about  rest  points  coinciding  with  local  minima  of 
the  averaged  potential  as  introduced  in  [7].  A  geometric  averaging  theory  of  mechanical  system  responses 
to  oscillatory  forcing  in  more  general  settings  (e.g.  capable  of  treating  the  “hovering”  motions  identified 
in  [7])  remains  to  be  fully  developed.  It  also  remains  to  explore  the  role  of  energy  methods  in  developing 
a  constructive  controllability  theory  along  the  lines  of  [10]. 
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